Correlation in narrow nanorods: a variational potential–configuration interaction scheme

Full configuration interaction calculations for two electrons in narrow semiconductor nanorods are carried out employing different orbital basis sets. It is shown that the usual configurations built from single-particle states cannot yield a correct singlet–triplet energetic order regardless of the basis size, as they miss the correlation energy. Mean-field optimized orbitals partially correct this drawback. A new approach is introduced, based on a simple variational procedure, which yields robust results.


Introduction
Quantum dots (QDs) are artificial semiconductor structures with sizes a few tens of nanometers. They are often modeled as local perturbations of the periodic crystal field of the semiconductor surrounding matrix, using envelope function and effective mass approximations to describe their electronic structure [1,2]. Within this approach the details of the unit cell are integrated out, so that only a macroscopic (or envelope) description of the system remains. Thus, the interaction between conduction electrons and atomic core electrons and nuclei is averaged as the interaction with a continuous polarizable medium. Likewise, Coulomb interaction between carriers is assumed to be screened by the dielectric constant of the QD bulk material [3]. The large value of the semiconductor dielectric constant coming into the Coulomb term (typically [10][11][12][13][14], along with the small value of the electron effective mass coming into the kinetic energy term, often lead to situations in which the separation between discrete single-particle (SP) levels exceeds by far the characteristic interaction energy between particles. As a consequence, SP descriptions of the many-body problem may be used, treating the Coulomb interaction between conduction electrons as a perturbation [4][5][6].
A better account of the two-particle interaction can be achieved through diagonalization of the many-body Hamiltonian in the basis set of non-interacting configurations. This is the so-called full configuration interaction (FCI) method, which plays a central role in quantum chemistry [7]. FCI provides benchmark results for the ground state energy and wavefunction, as well as for the description of excited states, as it is inferred from the McDonald theorem [8]. The drawback of the approach is its computational cost, originating from its often slow convergence. The key for convergence is the selection of a suitable monoelectronic basis set. As showed by Löwdin and Shull [9], the basis of natural orbitals is the one which requires the fewest configurations to achieve a given accuracy in the energy. Natural orbitals are the orbitals that diagonalize the one-particle density matrix. Therefore, exact natural orbitals should be extracted from the FCI wavefunction. Alternatively, one may use singleparticle (SP) orbitals, i.e. the eigenfunctions of the onebody operator coming into the many-body Hamiltonian. This approach is quite common in few-electron and excitonic QD calculations, where convergence problems do not usually arise [4,5,10].
In this paper, however, we identify a QD system where correlation energies are so strong that the basis choice becomes critical for overcoming slow convergence issues. Namely, we study narrow rod-shaped QDs (nanorods, NRs) [11,12]. As recently shown, the large aspect ratio of these semiconductor structures is responsible for some remarkably strong Coulomb interactions [13][14][15]. Here we demonstrate that, even in the simplest case of two interacting electrons, the CI method based on SP orbitals (SP-CI) fails to converge. Indeed, for long rods it predicts a triplet ground state, thus violating a fundamental theorem for two-body exact functions [16]. To solve this problem, we follow quantum chemistry CI methods and replace the regular SP orbitals by effective orbitals obtained through mean-field calculations, either Hartree-Fock (HF-CI) [7] or Kohn-Sham density functional theory (KS-CI) [17]. The suitability of using these mean-field optimized orbitals is discussed and compared with a new, simpler approach based on a two-fold variational procedure.

Theory
The Hamiltonian of interacting conduction electrons in a QD reads (1) where T = − 1 2m * ∇ 2 is the kinetic energy operator, v c the spatial confining potential, m * the effective mass and the dielectric constant.
In order to carry out a CI calculation, one first has to select a given one-electron basis set {φ p } M p=1 and then build all possible N-body Salter determinants α = det{φ 1 φ 2 · · · φ N } out of them. In a second step, the projection of Hamiltonian (1) onto the basis set of Slater determinants (or onto a subset of spin-and symmetry-adapted configurations) is carried out yielding matrix elements α |H| β . The last step is the diagonalization.
For practical purposes the one-electron basis set is taken as orthogonal, although many-body approaches based on overlapping orbitals, as for example the generalized valence bond (GVB), have been successfully developed in quantum chemistry [7,18]. A possible orbital choice is the SP states, i.e. the eigenfunctions of the one-electron Hamiltonian h. We must select a finite basis set. The natural choice is the M lowest-energy eigenvectors. This truncation of the basis set implies a projection of Hamiltonian (1) onto the corresponding FCI space yielding a model Hamiltonian that, in second quantization language, reads where ε p = h pp is the pth eigenvalue of h, g pqrs is the two-electron integral and a + p , a q the usual fermion creator/annihilator operators.
Alternatively, we may employ another M-dimensional orthogonal basis set {χ p } M p=1 . In this case, Hamiltonian (1) turns into where off-diagonal h pq terms arise because χ p is no longer an eigenfunction of h.
We may consider {χ p } as related to a given v m f (r ) potential operator which is added and subtracted to Hamiltonian (1), Now we determine {χ p } from the eigenvalue equation f (r )χ p (r ) =ε p χ p (r ) and rewrite equation (3) as HF, DFT and, in general, any mean-field calculation are particular cases of this procedure. For example, in HF, Within the framework of this procedure, our approach to deal with strongly correlated regimes is the use of a variational potential, i.e. a potential depending on parameters that allow optimization of the orbital basis set in the CI process.
The advantage of adding a mean-field potential is that, contrary to the SP orbitals which disregard electron-electron interaction, HF orbitals include it somehow averaged in the Coulomb and exchange terms, and KS orbitals further account for some of the correlation through the corresponding functional. This has proved useful in describing many-electron QDs [19,20]. In our approach, all terms in Hamiltonian (1) are treated on an equal footing to optimize the potential and hence the orbital basis set.
It is worth noting that the projected Hamiltonian (3) is not defined by the selected orbital basis set but rather by the linear space that it generates. Thus, for example, in [19], where HF spin orbitals of an N-electron system are expressed as a linear combination of a set of K SP states, a FCI expansion which employs either K -dimensional basis set yields identical result. The advantage of the strategy presented in [19] relies on a further truncation of the HF basis set.

Results and discussion
We carry out calculations for the lowest-lying singlet and triplet states of a two-electron semiconductor CdSe NR. The rod is composed of a cylinder of radius R = 20Å and variable length L, attached to two hemispherical caps of the same radius R at the extremes. The material parameters are effective mass m * = 0.13, dielectric constant = 9.2 and confining potential v c = 4 eV. The dominant electronic configurations for the studied singlet and triplet states are σ 2 1 and σ 1 σ 2 , respectively, where σ 1 , σ 2 are the lowest-lying and first excited orbitals with zero z-component of the angular momentum (M z = 0). We have checked that, as expected, already in the case of a short NR, many σ (M z = 0) SP states lie energetically below the first π (M z = ±1) energy level, and that higher angular momentum states are much more excited. Then we employ a basis set 1-n-1 in the CI calculations, i.e. a basis set including the n lowest-lying σ and the lowest-lying degenerate pair of π orbitals, the last orbital pair making, in all studied cases, a very small contribution to the singlet and triplet CI wavefunctions. Figure 1 shows a series of SP-CI calculations (basis set 1-6-1) of the singlet and triplet states of a R = 20Å two-electron semiconductor CdSe NR versus L. At first glance, the results seem reasonable. At L = 0 (a sphere) the singlet state lies energetically far below the triplet state, and the contributions to the CI expansion of configurations other than the dominant one are extremely small. As L increases, the energetic singlet-triplet gap TS decreases, the contribution of excited configurations to the CI wavefunction increases and, in the large L limit, a singlet-triplet degeneration occurs. However, a closer view reveals that beyond L = 170Å the triplet state becomes the ground state (see inset in figure 1). This fact is a violation of a known theorem, attributed to Wigner, stating that any general two-particle Hamiltonian including a symmetric and real potential energy operator has a singlet ground state [16]. Therefore, the result is qualitatively wrong. One may relate the failure to the truncation of the orbital basis set. In order to check it, we ran a series of SP-CI calculations of the singlet-triplet energy gap TS in an L = 200Å NR versus the size n of the (1-n-1) basis set. The results, shown in figure 2, do not allow any finite n yielding the theoretically expected singlet ground state.
We now turn our attention to the wavefunction. The inset in figure 2 shows a profile of the FCI singlet (solid line) and triplet (dashed line) charge densities along the NR vertical axis. We can observe that the triplet and, to a lesser extent, the singlet state develop a valley in the electronic density profile around the rod center. The presence of this valley allows the interacting electrons to reduce the repulsion energy, hence stabilizing the states. This groove is already present in the independent-particle description of the triplet state (through the σ 2 orbital), but not in the singlet ground state, where it is a pure correlation effect [15]. In the case of a poor correlation description, the singlet valley is not deep enough and the triplet is more stable, leading to the triplet ground state predicted above.
Since the singlet state has gerade symmetry and only double excitations to ungerade orbitals may contribute to the formation of the central valley in the density profile of the singlet, we explore the possibility of carrying out a CI expansion using gerade orbitals with an in-built valley. A first attempt is to follow a recent approach proposed by Abolfath and Hawrylak of using HF orbitals [19]. However, the approach is useless in our case, because in this scheme HF orbitals are expanded in terms of a (large) SP basis set and, as shown in figure 2, even an extremely large (finite) SP basis set cannot prevent the singlet-triplet reversal. Note that gerade/ungerade HF orbitals would be linear combinations of gerade/ungerade SP orbitals. Therefore, no gerade HF orbital developing a valley can be obtained as a linear combination of SP orbitals. We then take a different approach and carry out numerical HF calculations, build a many-electron Hamiltonian (5) out of the M lowest-lying HF eigenvectors and repeat the calculations shown in figure 1 with this new orbital basis set. In particular, the results collected in figure 3 correspond to a 1-6-1 HF basis set. One can see that HF-CI outperforms SP-CI. Indeed, at L = 170Å, HF-CI clearly predicts a singlet ground state. All the same, a singlet-triplet crossing is observed shortly after L = 200Å, indicating that the basis set is still insufficient.
We investigate in some more detail the large L region. To this end, we consider a L = 200Å NR, and carry out calculations employing a large basis set (1-10-1). In the first instance we run mean-field HF and KS (in the local density approximation) calculations. The corresponding meanfield potentials are represented in figure 4 left and center, respectively. We see that both have a similar shape, i.e. a profile versus z resembling a symmetric hill with the top at the rod center. In figure 4 (right) we show a simplified potential v tri , composed of a triangle of height H along the NR z-axis. All three potentials have a common feature: they destabilize the central region of the NR. In the case of HF and KS, the meanfield potentials are optimized in a self-consistent procedure. As far as the third potential is concerned, we optimized the height H of v tri variationally, i.e. we perform FCI calculations with states obtained as eigenvectors of a modified SP Hamiltonian, which adds v tri (H ) to the confining potential, and then select variationally the best H . This optimization is shown in figure 5 for a n = 10 basis set. We represent singlet and triplet FCI energies versus the triangle height H . For comparison, horizontal lines representing HF-FCI and KS-FCI are enclosed. Several relevant conclusions that may be drawn from this figure: (i) all three potentials recover the correct singlet/triplet energetic ordering; (ii) orbital optimization affects the singlet state far more than the triplet; (iii) from the variational principle and McDonald's theorem, we conclude that v KS outperforms v HF and the variational potential v tri outperforms both v HF and v KS . The good performance of v KS is remarkable in view of its typical failures in very-few-electron systems [21]. Also, we note that in this system HF-CI gives better results than the SP-CI scheme, contrary to the expectations of [19] for a small number of electrons. The singlet-triplet TS gap yielded by v HF , v KS and v tri through the FCI calculations are 0.2, 1.4 and 3.5 meV, respectively.
As pointed out in section 1, the correlation energies in the QD system under study are so strong that the choice of   figure 5. Solid, dashed and dotted lines are used for the lowest, first excited and second excited orbital, respectively. a single-particle basis set becomes critical for a correct CI description. This is nicely shown in figure 6, where the ρ = 0 z-profiles of the orbitals with largest contribution  to the CI are illustrated. The lowest-lying (σ ) SP orbital shows a prominent maximum at the NR center. Therefore, the a priori dominant σ 2 configuration of the lowest-lying twobody singlet state incorporates a strong Coulomb repulsion. Meanwhile, the same figure shows that HF and KS optimized orbitals present a central plateau instead. The observed plateau is the result of the mean-field optimization which, as pointed out in section 2, accounts for two-body interactions in an averaged way. In our approach, all terms in the many-body Hamiltonian are treated on an equal footing to optimize the orbital basis set and, as shown in figure 6, the result of the CI optimization yields a lowest-lying σ orbital with a valley at the NR center. Consequently, the σ 2 configuration already provides a qualitatively correct description of the central valley displayed by the two-body singlet ground state (see inset in figure 2). Similar arguments can be made for the low-lying excited orbitals, also shown in figure 6.
In the next step, in order to show the performance of our approach compared with standard mean-field calculations, we consider again the 1-6-1 basis set employed in figures 1 and 3, run calculations optimizing v tri in the region 100Å L 250Å, and compare the results yielded by v tri with the SP-CI and HF-CI methods. The combined results are summarized in figure 7, where the FCI singlet-triplet gap TS obtained from SP, v HF and v tri versus L, is represented. One can see that only in the case of v tri is the correct null asymptotic limit of the singlet-triplet gap obtained. In the inset, the difference between the calculated SP or HF gaps and the benchmark one (v tri ) is displayed. The inset shows that, within this region, the gap error is approximately constant for both approaches. It mainly comes from a deficient description of the singlet state (see figure 5).
To close this section we would like to note that the poorly correlated description of SP-CI is exclusive to NRs with a large aspect ratio. As a matter of fact, SP-CI calculations for a R = 50Å, L = 200Å (total length 300Å) CdSe NR employing a 1-6-1 basis set yields the correct singlet-triplet order, contrary to the R = 20Å NR we have analyzed in detail.

Conclusion
We have shown that the FCI correlation energy of two electrons in quasi-one-dimensional NRs calculated with the usual orbital basis set of numerical eigenfunctions of the one-electron operator h arising in the many-body Hamiltonian (SP-FCI scheme) leads to the wrong singlet-triplet energetic order. Irrespective of basis size, the correlation energy of the NR is very sensitive to the choice of orbitals. Thus, mean-field-CI approaches that somehow take into account an averaged electron-electron interaction in the orbital optimization can partially correct this drawback. A new CI scheme is proposed in which the mean-field potential is replaced by a simple triangular potential, with the triangle height as a variational parameter. The latter approach outperforms mean-field-CI schemes, yielding robust results.