All articles published by MDPI are made immediately available worldwide under an open access license. No special
permission is required to reuse all or part of the article published by MDPI, including figures and tables. For
articles published under an open access Creative Common CC BY license, any part of the article may be reused without
permission provided that the original article is clearly cited.
Feature Papers represent the most advanced research with significant potential for high impact in the field. Feature
Papers are submitted upon individual invitation or recommendation by the scientific editors and undergo peer review
prior to publication.
The Feature Paper can be either an original research article, a substantial novel research study that often involves
several techniques or approaches, or a comprehensive review paper with concise and precise updates on the latest
progress in the field that systematically reviews the most exciting advances in scientific literature. This type of
paper provides an outlook on future directions of research or possible applications.
Editor’s Choice articles are based on recommendations by the scientific editors of MDPI journals from around the world.
Editors select a small number of articles recently published in the journal that they believe will be particularly
interesting to authors, or important in this field. The aim is to provide a snapshot of some of the most exciting work
published in the various research areas of the journal.
A multireference variant of coupled cluster theory is described that applies to systems that can qualitatively be described by deleting two electrons from a closed shell determinant, for example biradicals, single bond breaking processes, or valence excited states. The theory can be generalized to arbitrary open-shell systems and takes a form that is akin to equation-of-motion coupled cluster theory, but where all wave function parameters are explicitly optimized for the state of interest. The implementation of the present methods was accomplished in an automated fashion using the recently developed Automatic Program Generator (APG). We present benchmark results for the O2 and F2 molecules and investigate the behaviour of a number of closely related variants within the same general framework.
Single reference coupled cluster theory, [1,2] especially CCSD(T) [3,4], is a highly successful approach to describe the electronic structure of systems that are qualitatively reasonably well described by a single determinant (many closed-shell molecules and high-spin open-shell systems) [5,6,7,8]. Likewise coupled cluster response theory (CCLRT [9,10,11]) and its cousin equation-of-motion coupled cluster theory (EOM-CC [5,12]) are well suited to describe singly excited states and radicals. An important prerequisite for the applicability of the CC response approach is the presence of a nearby state that can serve as the reference in the parent state CCSD calculation. Most commonly the reference state is a closed-shell system that differs by up to one electron from the states of interest. The standard CCSD equations are solved for the reference state and the primary role of the cluster operator is to include dynamical correlation, which tends to be transferable from one state to the other. The electronic states of interest are obtained by diagonalizing the transformed Hamiltonian over the proper configurations. In practice, if an excited state is dominated by 1p1h (i.e. single) excitations (or 1h, or 1p substitutions in case of IP-EOMCC and EA-EOMCC respectively [13,14,15,16]) one needs to include one more excitation level in the diagonalization space to achieve an accuracy of about 0.1 - 0.3 eV [17,18]. The accuracy of EOMCC and/or CCLRT can be further improved by including perturbative triple corrections (of the iterative or non-iterative type) [19,20,21,22,23], while also full EOM-CCSDT results have been presented [24,25]. Recently, Krylov [26,27] pushed the EOM strategy a level further by choosing as the reference state a high spin (usually triplet) state and diagonalizing over a manifold of states of different spin multiplicity (the so-called spin-flip EOM method). If the dominant part of the final state wave function resides in the principal sector (related to the reference by a single spin-flip excitation) the approach has been found to yield promising results. The (triplet based) spin flip EOM approach provides access to biradicals and the breaking of single bonds and appears a nice extension of the tools in the quantum chemistry community.
Other CC variants exist to describe excited states that involve yet another similarity transform of the Hamiltonian related to ionization potentials and electron affinities. The purpose of the second transformation is to remove (or decrease) the coupling of the primary sector to the remaining sectors [28,29,30,31]. Fock Space Coupled Cluster [32,33,34,35] and Similarity Transformed EOMCC (STEOM [30,36,37,38]) are the best known variants of such an approach and they provide a very efficient approach to excitation spectra as the final diagonalization space can be very compact. Moreover, these approaches have traditionally also been used to describe certain multireference situations. A recent example is the description of the electronic states of the NO3 cation , most of which have significant multireference character. The NO3- anion is well described by CCSD and this is the state that is used as a reference, while the final diagonalization is over states that remove two electrons from the reference, e.g. . For this reason the method is referred to as the double-ionization potential STEOM (DIP-STEOM) approach. Beyond computational efficiency a strong feature of the FSCC and STEOM approaches is the fact that they are rigorously spin-adapted and size-consistent. A weakness is that the reference state may be far removed from the actual states of interest and orbital relaxation effects may be substantial which are difficult to describe through configuration mixing. The reference state can become entirely artificial in these schemes as exemplified by results for the vibrational frequencies of ozone using the DIP-STEOM approach based on the di-anion of O3
. In reality the dianion of ozone is non-existent and the bound reference state is an artifact of the basis set. Therefore in such cases the DIP approach can only give reasonable (and even excellent) results in small basis sets, lacking diffuse character in particular.
The idea behind all of these approaches is that information on dynamical correlation coefficients are obtained from states that are well described by single reference theories, and these are then transferred to the actual states of interest (through similarity transformations [28,29,30,31]). The general success of these approaches indicates that dynamical correlation is largely transferable (the concept of valence universality, e.g. [34,40]), and it also tells us that the parameterizations for the wave functions can be quite compact. All of the above approaches can be viewed as internally contracted methodologies [41}: The cluster operator acts upon the qualitatively correct wave function as a whole. However, the fact that in practice the cluster coefficients have to be obtained from a well-behaved single reference state is a limitation. Such a state is not always present, or its presence may depend on the nuclear geometry considered, while the transferability of dynamical correlation will depend on the similarity of reference and target states. These statements are equivalent to the fact that valence universal (or Fock Space) theories as discussed in detail by Mukherjee and Pal , and Paldus and Jeziorski  have their limitations for the general open-shell problem. The parameterizations are very suitable in my mind, but the equations to determine the amplitudes, i.e. the valence universality concept may break down.
In this paper we investigate simple multireference situations that can be characterized by the fact that one virtual orbital dominates the correlation effects. Two simple diatomic molecules, O2 and F2, are used as examples, and they share the property that qualitatively correct wave functions (for ground and valence excited states) can be described by creating two holes in a formal closed-shell di-anion determinant. Therefore our basic formulation follows a double ionization type of strategy. We will use the same parameterization for the wave function as is used in the double ionization variant of equation-of-motion coupled cluster theory (DIP-EOMCC ). The equations that determine these parameters are quite different however. In principle all parameters are determined specifically for the state of interest. We will discuss a number of slightly different variants and generalize the state-specific scheme such that coefficients can be determined in a state-averaged way. Finally we come full circle and investigate if the state specific coefficients for one particular state could also be used for other states, implying that they would again be transferable or valence universal. The approaches we discuss in this paper are steps towards a completely general internally contracted multireference coupled cluster approach along the lines discussed in [29,42,43]. However, none of the approaches we discuss here are rigorously size-consistent, in the sense that an open-shell system does not separate quite correctly into open-shell subsystems. For this reason we term our approaches state-selective variants of equation-of-motion coupled cluster theory, rather than multireference CC. The methods are size-intensive in the sense that the energies scale properly if closed-shell systems are added at infinity.
In this paper we consider wave functions that are qualitatively well described by the parameterization
The closed-shell determinant contains two more electrons than the state(s) of interest and serves as the vacuum in the many-body theory. Occupied orbitals in are indicated as i,j,k,l etc., while virtual orbitals carry labels a,b,c,d. The state is the reference state and is assumed to contain the dominant configurations in the target wave function. The orbitals in are defined such that the reference state is optimal in some sense. The orbitals can be defined for example using an MCSCF procedure. Let us note here that our reference state has a slightly different character than in other multireference approaches, as there are no active space restrictions in Eqn. (1). Rather, the 'active' space in our formulation is implicit by the definition of the vacuum state that contains one more spatial orbital than the number of electron pairs. Once the closed-shell vacuum determinant is defined all orbitals are treated equivalently.
The complete correlated wave function is parameterized as
is the traditional closed-shell CCSD operator, where particles and holes are defined with respect to . The linear operator has the form
such that a complete operator (in a full CI sense containing up to (N + 2)h − Np excitations) can describe any state with 2 less electrons than . Substituting the parameterization in the Schrödinger equation and multiplying by we find
The operator is hence the eigenvector of a transformed Hamiltonian, as in EOMCC. In principle the transformation with can be arbitrary, and its purpose is to improve the convergence properties of the truncated diagonalization, that we wish to carry only to the 3h1p level. The arbitrariness or complete redundancy of the operator in the exact case (i.e. full CI expansion of ) indicates that a proper definition of the equations for is not straightforward. It also implies that different values for could lead to identical forms for the wave functions if the operator is adjusted accordingly. This arbitrariness/redundancy is at the heart of the difficulty of formulating a multireference CC theory. Formal theoretical considerations can help only so much. Actual computational schemes have to be tried out in practice and their robustness and general applicability have to be established by elaborate testing.
While the proper equations for require further discussion (see below) the equations for the operator , truncated as in Eqn. (4) are clear-cut. We will simply diagonalize over the configuration space defined by , hence
In this work the equations for are defined in a state dependent way by projecting the Schrödinger equation (5) against excitations out of the reference state, , where . In practice we usually define the reference state through , but the theory will be developed more generally, and can in principal be obtained from a separate calculation, of MCSCF type for example. The projected equations read
The analogous equation for would be entirely redundant however, as is entirely contained in the manifold of Eqn. (6). After some experimentation we found that the following is a suitable alternative
in which we essentially replaced by the dominant part, the reference state . The resulting equations can still exhibit exact redundancies or near singularities. This implies that quite different values for the parameters can generate (nearly) the same wave function, and it typically means that the non-linear equations are very hard to solve numerically. Similar issues arise in internally contracted MR-CI  and in practice we eliminate the near redundancies, by projection on a regularized subspace. This is non-trivial and requires further discussion (see below). The above equations are coupled and we solve for the , amplitudes from one iterative procedure. This means that each state has its associated set of amplitudes and the procedure is completely state-specific. In addition we have the flexibility to use the amplitudes (and t1-equation) to define Brueckner orbitals.
We will also consider an alternative and simplified scheme where we define the t2-equations as
in analogy to the equation for the t1 coefficients.
The equations in which we use the complete operator in the t2-equations (Eqn. 7) define the Complete State Selective EOMCC approach (CSS-EOMCC), while we use the acronym RSS-EOMCC if instead of we use (Eqn. 9). As a mnemonic the acronyms contain the respective operators or , while in addition the abreviation RSS-EOMCC stands for "Restricted (or Reference) State Selective EOMCC". While in the CSS-EOMCC approach , in the RSS-EOMCC approach we take an additional short-cut by obtaining from a diagonalization of the bare Hamiltonian over the 2h configurations only (DIP-TDA), and hence is not relaxed in the presence of dynamical correlation effects. Both of these methods are initial trial versions and in this paper we simply want to investigate the sensitivity of the results to such minor variations on the theme. There are a number of details that require further discussion. Let me briefly mention the salient points and defer the full details to a future exposition. In due time I hope to establish a single robust approach in which these details are specified once and for all. It does not serve much of a purpose to invent a large set of slightly different variants except in the testing stage.
In the above formulation we used spin-orbitals to facilitate the discussion somewhat but this is not what is done in practice. In reality we parameterize the operator using generators of the unitary group
and the projections used in the T-equations are against and . The parameterization of the operators in terms of the generators of the unitary group ensures that the transformed Hamiltonian commutes with the spin-operators and the resulting wave functions are rigorous spin eigenfunctions. Let us note that the spin-adapted variant of the equations are not equivalent to the spin-orbital equations in this case. These issues have been discussed in detail in a previous paper of us in which we explore spin-adaptation in a CC framework . The final diagonalization of over the 2h and 3h1p configurations is carried out in a spin orbital basis but the diagonalization manifold is easily seen to be spin adapted. In principle we could use an explicitly spin-adapted formulation in this step, and spin-orbitals are used purely for our convenience.
Elimination of near-singularities.
The reference state defines a symmetric and positive definite metric matrix and the doubles analog , neither of which depend on the virtual labels a,b. These metric matrices have orthonormal eigenvectors Uλ and eigenvalues λ, some of which may be small (smaller than 0.01 say). Defining projectors on the regularized subspaces we solve the regularized equations
with similar expressions for the singles equations. The reference state entering the projectors is determined at the beginning of the iterative procedure and is kept fixed. This improves convergence, compared to projections that are determined iteratively. We have an input threshold variable to decide which eigenvalues of the metric matrix are kept in the projection, and this is set to 0.01 in the current calculations.
Size consistency, size-extensivity and size-intensivity.
If a compound system consists of an open-shell system A and an additional non-interacting closed-shell molecule B (located at infinity), there exist solutions to all of the equations in this paper that correspond exactly to the subsystems in isolation. In particular, if the molecular orbitals are taken to be localized all operators , , become completely localized on the individual subsystems. The localized equations for closed-shell subsystem B in the compound system reduce to the standard CCSD equations for system B alone. Likewise, the localized equations for open-shell subsystem A are unaffected by the presence of non-interacting subsystem B. In this case of a compound system consisting of an open-shell and a closed-shell system the orbitals can always be localized and this facilitates the analysis. However the theory is invariant under such localization and the conclusion is completely general. This is a very desirable scalability property that we refer to as size-intensivity. It presumably implies that the energies, in particular energy differences between different open-shell states, are size-extensive (= physically well-behaved) as the size of the system grows, although a precise mathematical definition of this property for an open-shell system is hard to give. Let me emphasize this point here. To my knowledge there is no theoretical framework to analyze or even define the property of size-extensivity in the thermodynamical sense for an open-shell system, because there is no repeating unit that leads to a constant energy per unit cell in the infinite limit.
While the SS-EOMCC approach separates properly into an open-shell plus closed-shell fragements (size-intensivity) it is certainly not size-consistent in the sense that if a triplet molecule separates into two doublet radicals, the energy we would calculate using the present method is not precisely twice the energy of the doublet states in isolation (calculated by the analogous method for doublet states ). This is due to the linear CI diagonalization for the target states. The origin of the problem is the same as the issue with charge-transfer separability in EE-EOMCC calculations as analyzed in detail in earlier papers [30,44,45].
The details of the equations (in particular the equations for the t2 manifold) are quite involved. When we started this work we were quite aware that we needed to be able to explore various alternatives for a possible MRCC theory. Following earlier work in this direction, in particular by Janssen and Schaefer  and Li and Paldus  we built a tool to aid the derivation and implementation of many-body methodology using automation. The Automatic Program Generator (APG) we developed, generates an efficient Fortran procedure to perform the calculation starting from simple input equations like the ones used in this paper. Some details can be found in our earlier paper on a MRCC theory for doublet states . The APG derives the detailed equations using essentially Wick's theorem and some further embellishments (e.g. taking symbolic derivatives, multiplying by density matrices, spin-integration). The resulting equations consist of a large number of terms (on the order of 300 spin-adapted Goldstone diagrams), consisting of multiplications of typically 3 or 4 operators. In the next step the equations are factorized so that each elementary step in the factual calculation consists of addition or contraction of two matrices. Optimal factorization constitutes an n-p complete problem in mathematics, and at present we appear to have developed a reasonable strategy to tackle the factorization issue. Finally the APG generates a Fortran subroutine that calculates the residuals of Eqns. (7-9) given a set of input amplitudes, the relevant density matrices of the reference state, and Hamiltonian matrix elements. The generated Fortran subroutines use a library of handwritten subroutines that are computationally efficient. This efficiency is illustrated by the fact that the same subroutines have been used in STEOM-CCSD calculations on free base porphyrin in a DZP+diffuse basis set . The APG is well tested [49,50] and can be expected to provide fairly efficient and error-free implementations of highly advanced electronic structure methods.
An additional factorization step.
In order to facilitate the computation, in particular in connection to the APG we introduce an additional approximation in the case of a CSS-EOM-CCSD calculation. In practice we factorize the operator as , where , where m indicates an active occupied orbital label. The operator has the same form as the IP-operator in STEOM-CCSD theory (e.g. ). This factorized form of only occurs in the evaluation of the t2 equation (Eqn. (7)), not in the final diagonalization. The reason to use this factorization step is that we can then avoid the use of three-particle intermediates in any of the factorized equations, and we can use unitary group operators throughout. Let me emphasize here that this step has a purely computational origin, and it is essential to employ this factorization in conjunction with the present version of the APG. We do introduce an active space in this step, but this is only done to fit the true operator. The factorization is a computational device, rather than a theoretical construction, and if the CSS-EOMCC scheme were to turn out our final method of choice we can eliminate this factorization step with some additional work. In practice we solve a regularized form of the following linear equation
to obtain the coefficients of the operator. From our results below we find that the contribution due to itself only marginally affects the t2 amplitudes, and hence the slight error we make in the above factorization is presumably of little importance. The choice of active space is a purely computational matter, and we simple choose it to optimize the fit. It has nothing to do with the multiconfigurational problem.
Truncation of singles.
In practice we do not fully develop the quartic expansion of the single excitation coefficients, but instead only keep up to linear terms. The reason for this is again purely technical. If Brueckner orbitals are used the results are unaffected by the truncation. For the examples considered in this work the truncation makes very little difference, as the t1-amplitudes are always small (0.02 or less).
State averaged Calculations.
It can be advantageous to perform calculations in a state averaged way. This is true in particular if degeneracies are involved (e.g. Π and Δ states for the diatomics below). The state selective determination of the t-amplitudes can break symmetry, so that diagonalization of would not lead to equal energies for the Πx and Πy states in a diatomic for example. In a state averaged calculation only the equations for the t-amplitudes are modified:
In principal we have a flexibility to choose the weights wλ in the sum over included states, although in practice we always keep them equal (essential to keep degeneracies). The implementation of the programs is flexible and the subroutines written by the APG take the state-averaged density matrices on input.
All of the implementations were incorporated in a local version of the ACES II package , which does not have any MCSCF capabilities at present. We used the APG to implement a limited MCSCF procedure through a CC singles variant of our SS-EOMCC equations. In the Brueckner MCSCF scheme we solve the following set of coupled equations
where the diagonalization of is over the reference 2h configurations. The t1 amplitudes are used to define a rotation of the orbitals defining and the usual Brueckner scheme is used to iterate to convergence .
Iterative Sequence of the Solution Algorithm
In the first step of the algorithm we find a suitable set of starting orbitals to define . Typically we use UHF or ROHF orbitals for the triplet state, from which we construct a spatial orbital density matrix. The strongly occupied natural orbitals of this density matrix define the closed-shell reference state. This scheme is analogous to Pulay's UNO-CAS scheme , and the orbitals can be further refined using the Brueckner-MCSCF procedure. We can also start from converged RHF orbitals for the dianion. After the determination of the vacuum state the algorithm has two branches depending on whether we use the CSS-EOMCC or RSS-EOMCC variant.
Steps in RSS-EOMCC calculation:
Calculate , .
Solve the RSS-EOMCC equations for t1 and t2 amplitudes, using DIIS  to accelerate convergence.
Calculate appropriate matrix elements of (in particular the numerous and elements are not needed).
Diagonalize over the 2h and 3h1p configurations. We can easily find multiple states in this scheme. The operator is defined in a state selective fashion, but if the t-amplitudes are assumed to be transferable, a number of states can be calculated simultaneously. This latter scheme is referred to as t(ransfer)-RSS-EOMCC.
Steps in a CSS-EOMCC calculation
Calculate , .
Iterate: [ Calculate appropriate elements of ; Diagonalize over 2h and 3h1p configurations to obtain ; Fit to ; Calculate residuals in the and equations and update the t-amplitudes, using ]. In our DIIS extrapolation we use both the and operators (and orbital rotation parameters in a Brueckner calculation following reference ).
In a final diagonalization of we can include more states and in this way gain insight into the transferability of t-amplitudes. Such calculations are acronymed t-CSS-EOMCC.
The amount of work that goes into a CSS-EOMCC calculation is substantially greater than in a RSS-EOMCC calculation: In each CSS-EOMCC iteration we form and diagonalize over the 2h/3h1p determinants. Moreover the CSS-EOMCC t2-equation is far more involved due to the presence of (in the guise of ).
We will present results for ground and valence excited states for the diatomics O2 and F2. In both cases we used the cc-pVTZ basis set , and we compare a variety of coupled cluster based methods.
III.a Results for O2.
In table I we report the optimized internuclear distance and the harmonic vibrational frequency for the three valence states of di-oxygen, while in table II the 0-0 transitions between ground and excited states are reported (using the harmonic approximation). All three electronic states have essentially two electrons in the four πg spin-orbitals, and the mixing between the 6 determinants leads to 3 electronic multiplets , 1 Δg, . Only the triplet ground state can be accessed by conventional single reference methods and the CCSD(T) results in the cc-pVTZ basis set is seen to compare quite satisfactorily to experiment  (see table I). We have also included DIP-STEOM and DIP-EOM results, and we find that in particular the DIP-STEOM results are very satisfactory. In the DIP calculations we use the orbitals of the di-anion to build the reference state. This works satisfactory in the cc-pVTZ basis but results would deteriorate sharply if more diffuse orbitals are included in the basis set. In DIP-STEOM calculations "triple excitations" (4h2p type of configurations) are included implicitly and this presumably explains the slight edge that DIP-STEOM has over DIP-EOM. In table I & table II we included four variants of State Selective EOMCC results. In the t-RSS and t-CSS results we optimize the t-amplitudes for the ground state triplet and then use these amplitudes to calculate all three excited states of interest through diagonalization of . The vacuum determinant is built from triplet ROHF orbitals here. In all of the reported SS-EOMCC calculations the 1Δg state is obtained from a state averaged calculation that includes both partners of the multiplet. This avoids symmetry breaking. Analyzing the results, we find that the difference between the RSS and t-RSS results is completely negligible, and this indicates that the t-amplitudes are to a large extent transferable between these states that are similar in character of course. This conclusion does not hold very well, however, for the CSS and t-CSS results which exhibit significant differences. From a formal point of view the CSS results are always preferable over t-CSS of course, as all amplitudes in CSS are optimized for the state of interest. The difference between t-CSS and CSS is a little surprising however and it may signal that something is not quite stable in the CSS scheme. For the triplet ground state the only difference is the use of Brueckner orbitals in the RSS and CSS schemes, vs. ROHF orbitals in the t-RSS and t-SCC schemes. This does not make a significant difference. As far as geometries and frequencies are concerned the DIP, CSS, RSS and t-RSS methods perform about equally well, while t-CSS-EOM-CC is relatively worse.
For comparison we also included results from Brillouin-Wigner MRCC calculations, reported recently . The MR-BWCC results are not quite as good as the present results and we think this is likely due to the reference space used in the MR-BWCC calculations that only included the qualitatively essential determinants that differ in the occupation of the πg spin-orbitals. There is a substantial admixing from a doubly excited configuration in all of these states which is included in the reference state and the projection manifold in our DIP and SS based calculations, but not in MR-BWCC. We think that the choice of reference space is the dominant reason for the improvement of the present results over MR-BWCC, and the results probably do not reflect intrinsic differences in methodology. A further comparison is made with results from MR-AQCC calculations  in the cc-pVTZ basis set, where the 1s core orbitals were dropped from the calculation. The active space in the MR-AQCC treatment comprised the full valence space. The results for the triplet can be compared to dropped core ROHF-CCSD(T) and even ROHF-CCSDT calculations and we find that the MR-AQCC equilibrium bond distance deviates by about 0.5 pm. The effect of dropping the core is about 0.35 pm, and we find that for the triplet state the bond distance in the dropped-core MR-AQCC deviates by about 0.9 pm from the full CCSD(T) result. As seen in the table the difference between CCSD(T) and CCSDT is negligible, so we think the ROHF-CCSD(T) is a suitable reference point. For the other states the difference between our SS-EOM results and the MR-AQCC result is also around 1 pm, and so the differences in bond length correlate rather well. An analogous result holds for the vibrational frequencies. The MR-AQCC results are consistently about 80 cm-1 lower than the RSS-EOM-CCSD results, which appears rather much. The difference between MR-AQCC and CCSD(T) in the triplet state amounts to 35 cm-1, the difference due to the dropping of the core is about 15 cm-1, while on the other hand the RSS-EOM-CCSD results overshoot the CCSD(T) result by about 30 cm-1. These deviations all point in the same direction and lead to the rather large total deviation of about 80 cm-1. It is satisfying that this deviation is quite constant between the RSS-EOM-CCSD and the MR-AQCC approaches. The deviation appears slightly less systematic when comparing to the CSS-EOM-CCSD approach.
Comparing the calculated 0-0 transition energies (based on harmonic vibrational frequencies) in table II to experiment we find that the CSS-EOMCC results are excellent, differing by at most 0.015 eV. Results are overall pretty good however, and the errors do not exceed 0.1 eV. MR-AQCC results, extracted from reference  have a similar quality, while the MRBW-CCSD results are somewhat worse for the transition into the state, where the error amounts to 0.28 eV.
III b. Results for F2.
In our second set of test results we consider the ground and valence excited states of the fluorine molecule. In table III we have collected ground state properties for various methods. The CCSD(T) method serves as a benchmark and a comparison with experiment  shows that the basis set is far from being saturated. Comparing our results therefore with CCSD(T) rather than experiment we find that both of the DIP methods, which are based on an artificial di-anion reference state, yield rather poor results. Moreover, these results can easily be made far worse by selecting a basis set with diffuse functions. The results from the CSS-EOMCC calculations on the other hand compare well with CCSD(T), and we note that the RSS-EOMCC results are even slightly better. The total energy is significantly better (by 7 mH) in RSS-EOMCC compared to CSS-EOMCC.
In the DIP based scheme we can access single and double excitations of fluorine into the low-lying σu LUMO orbital. In this case we cannot compare to experimental results as all of these excited states are dissociative, while also the basis set is not particularly suitable for excited state calculations. Some of the triplet excitation energies can be obtained from ΔCCSDT calculations and they can serve as a benchmark. All calculations have been performed at the CCSD(T)/cc-pVTZ ground state geometry (R=141.36 pm). The singlet and triplet states of a given symmetry have been obtained from state averaged calculations, while in addition we always included complete multiplets in the state-averaged calculations in case of spatial degeneracies (Π and Δ states).
We have collected the orbital characters of the various singlet and triplet states in table IV. Here we indicate the orbitals that are annihilated from the fully filled shell defining F22- which is isoelectronic with Ne2. States that do not contain a deletion involving the σu LUMO are in fact doubly excited with respect to the ground state (i.e. 1Δg, 3Δu, , , ). Most states (including the ground state) have a fair amount of double excitation character. The singly excited states can be calculated by conventional EE-EOM-CCSD and EE-STEOM-CCSD methods, while the DIP and SS based methods can describe all valence excited states. From tables V (singlet states) and table VI (triplet states) it is seen that quite consistent results are obtained. The agreement between RSS-EOM-CCSD and CSS-EOM-CCSD is good, typically better than 0.1 eV. Moreover for the single determinantal excited triplet states (3Πu, 3Πg, ) we can compare to ΔCCSD(T) results and also here the comparison with SS-EOMCC is excellent. Since there is little reason to expect some states to be significantly more accurate than others in the SS-EOMCC methods we think that the accuracy is probably good across the board, although we do not have benchmark results to substantiate this assertion. The results from the DIP-STEOM approach are consistently low compared to the other schemes. This might reflect implicit triple excitation effects that are only present in DIP-STEOM, but given the rather poor comparison to CCSD(T) (which includes triple effects explicitly), it appears that the ground state is perhaps somewhat high in energy compared to the excited states in DIP-STEOM. Somewhat surprisingly the EE-EOMCC and EE-STEOMCC methods provide results that can be in significant error (from about 0.2 eV up to 0.5 eV for the state). This is presumably related to the fact that the ground state of F2 is fairly highly correlated (largest t2 amplitude is -0.18): large electron correlation effects in the parent state in general tend to deteriorate results in EOMCC and STEOMCC calculations. The mediocre results from the not-state-selective EOMCC and STEOM calculations indicate the difficulty to obtain a balanced description of ground and excited states of the fluorine molecule, and provide some perspective as to the high quality of the SS-EOM results. Finally let us discuss if the t coefficients are transferable from one state to the other in an SS-EOMCC calculation. Results are presented in the last three rows of tables V and VI. We have taken the t coefficients of either the state or the ground state t-coefficients and diagonalized the corresponding transformed Hamiltonian. Very clearly the results are not good at all in this case (errors of 0.5 to 1.5 eV occur). This is a little surprising as these kind of effects do not seem to occur in regular equation of motion coupled cluster calculations or in spin flip EOM as investigated recently by Krylov . One possible explanation is that some correlation effects may be missing in the transfer scheme. For example in the calculations based on the reference state there are no double excitations that involve the σg orbital as it is half empty, while such excitations should be present in the Π and Δ states, but using the cluster amplitudes based on they are not. Another factor may be that underlying the standard EOM schemes is a many-body similarity transform. Hence the contributions in the final diagonalization from more highly (and neglected) excited determinants can rather easily be analyzed (e.g. [31,49]) through the couplings present in the transformed Hamiltonian. In the current state-selective scheme this strict many-body analysis is discarded and this may well be the source for the lack of transferability of dynamical correlation effects. I regard this aspect a little suspect and hope to investigate the issue more deeply in the future.
There is one more topic that warrants some further consideration: total energies. In Table VII we have listed the partitionings of the total energies for the CSS-EOM-CCSD and RSS-EOM-CCSD calculations. We included the states for O2 and the singlet states for F2, (the data for the triplet states of F2 is similar due to the use of the state-averaging strategy). The energy is partioned according to the vacuum determinantal energy corresponding to the formal dianion, the closed shell part of the correlation energy (in a diagrammatic sense) , and the total energy. For comparison in table VII we listed results for Brueckner based CSS-EOM-CCSD and RSS-EOM-CCSD evaluated at identical geometries. First of all we note that even though the formal reference determinant is the same for all states, there are significant differences in vacuum energy that are due to orbital relaxation effects. So even though all states in dioxygen have the same spatial configuration, this does not mean that the optimum spatial orbitals are the same for all states. The ground states in both O2 and F2 have the highest vacuum state energy (corresponding to the dianion), and this reflects the relative compactness of the ground state LUMO. The closed-shell parts of the correlation energy tend to be fairly similar, a clear exception being the closed-shell part of the correlation energy of the ground state for F2. There are significant differences however between the CSS and RSS variants, in particular the difference in vacuum state energy corresponding to the respective Brueckner determinants. The RSS-EOM-CCSD method universally yields lower total energies than CSS-EOM-CCSD, on the order of 10 mHartree (0.3 eV). The relative energies between various excited states are far more accurate than 0.3 eV, but it does remain a little puzzling where the discrepancy in total energy comes from. Part of it seems to be related to the vacuum energy part which can easily be different by 25 mHartree. It may be good to reiterate that in the RSS approach the reference state is obtained from a DIP-TDA calculation, diagonalizing the bare Hamiltonian over the 2-hole configurations, while subsequently the t-amplitudes are solved and the vacuum determinant is defined in a Brueckner sense. In contrast in the CSS-EOMCC approach at each step the DIP-EOMCC equations are solved and the reference state is defined as the 2-hole component of the full DIP-EOM eigenvector. So in the CSS-EOMCC approach the reference state is fully relaxed with respect to dynamical correlation. It is counterintuitive that the RSS-EOMCC approach would yield the lowest energy and this is certainly a clear manifestation of the lack of a variational principle associated with these approaches. The comparison also makes clear that to obtain accurate energy differences we are relying on a balance that leads to a cancellation of errors. I have little insight in the systematics of this balancing act at present and would certainly have preferred to get more consistent total energies.
In this paper we have presented results for a few state selective variants of equation of motion coupled cluster theory (SS-EOMCC). In principle all parameters: orbitals, "closed-shell" cluster coefficients and wave function amplitudes of the transformed Hamiltonian, can be optimized explicitly for the state of interest. We limited the discussion here to situations were the dominant part of the wave function can be described by deleting two electrons from a closed-shell vacuum determinant, the so-called DIP scheme. The scheme can be readily generalized to the deletion of n-electrons from an (N+n)-electron type closed-shell vacuum state, and this appears to be a promising step towards a completely general internally contracted multireference CC method. For such a general scheme it is clear that a state selective scheme is required as the vacuum state loses all physical significance. We discussed two principal variants of the DIP based SS-EOMCC scheme. In the Complete SS-EOMCC scheme the t2 amplitudes are obtained by projecting the parameterized Schrödinger equation onto a suitable excitation manifold, while in the Reference SS-EOMCC scheme we use only the dominant 2-hole part of state and project , where is obtained from a DIP-TDA calculation. The latter scheme is quite a bit simpler and it provides equally good, possibly better results than the complete variant, in particular for the absolute total energy. Both methods overall perform quite satisfactorily for the diatomic molecules considered, especially for energy differences, and it is premature to make a definitive assessment or comparison between the two approaches. More results and benchmarks need to be obtained. However, if the difference between the two approaches is consistently minor the RCC-EOMCC method is definitely to be preferred, because of its relative simplicity. Let us emphasize here that at present we do not have a rigorous formal argument to choose one approach over the other, as the amplitudes are abitrary in the limit of a complete expansion. We also investigated if the t-amplitudes that are obtained in this scheme are transferable from one state to the next, and this appears not to be the case. This is a slightly surprising result as the t-amplitudes are supposed to account purely for dynamical correlation, which is expected to be fairly similar between related states. The observations made in this study may guide us in further explorations in our quest for completely general multireference CC methods.
I gratefully acknowledge NSF funding for this work through the Information Technology Program. I would also like to thank Prof. Petr Čárski for sharing the MR-BWCC results on singlet oxygen prior to final publication.
Coester, F. Nucl. Phys.1958, 7, 421.
Cizek, J. J. Chem. Phys.1966, 45, 4256.
Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. Chem. Phys. Lett.1989, 157, 479.
Bartlett, R. J.; Watts, J. D.; Kucharski, S. A.; Noga, J. Chem. Phys. Lett.1990, 165, 513.
Bartlett, R. J.; Stanton, J. F. Rev. Comp. Chem.1994, 5, 65.
Bartlett, R. J. Advanced Series in Physical Chemistry; Yarkony, D. R., Ed.; World Scientific, 1995; Vol. 2. [Google Scholar]
Lee, T. J.; Scuseria, G. E. Quantum Mechanical Electronic Structure Calculations with Chemical Accuracy; Langhoff, S. R., Ed.; Kluwer Academic Publishers: Dordrecht, 1995; pp. 47–108. [Google Scholar]
Mukherjee, D.; Pal, S. Adv. Quantum Chem.1989, 20, 292.
Kaldor, U. Theor. Chim. Acta.1991, 80, 427.
Nooijen, M.; Bartlett, R. J. J. Chem. Phys.1997, 106, 6441–6448.
Nooijen, M.; Bartlett, R. J. J. Chem. Phys.1997, 106, 6449–6455.
Nooijen, M. J. Phys. Chem. A2000, 104, 4553–4561.
Wladyslawski, M.; Nooijen, M. ACS conference proceedings, in press2001.
Jeziorski, B.; Paldus, J. J. Chem. Phys.1989, 90, 2714.
Knowles, P. J.; Werner, H. J. Theor. Chim. Acta1992, 84, 95–103.
Nooijen, M.; Bartlett, R. J. J. Chem. Phys.1996, 104, 2652–2668.
Nooijen, M.; Lotrich, V. J Mol. Struc.-Theochem2001, 547, 253–267.
Koch, H.; Jensen, H. A.; Jørgensen, P.; Helgaker, T. J. Chem. Phys.1990, 93, 3345.
Stanton, J. F. J. Chem. Phys.1994, 101, 8928.
Janssen, C. L.; Schaefer, H. F. Theor. Chim. Acta1991, 79, 1–42.
Li, X.; Paldus, J. J. Chem. Phys.1994, 101, 8812.
Gwaltney, S. R.; Bartlett, R. J. J. Chem. Phys.1998, 108, 6790–6798.
Nooijen, M.; Lotrich, V. J. Chem. Phys.2000, 113, 494–505.
Nooijen, M.; Lotrich, V. J. Chem. Phys.2000, 113, 4549–4557.
Stanton, J. F.; Gauss, J.; Watts, J. D.; Nooijen, M.; Oliphant, N.; Perera, S. A.; Szalay, P. G.; Lauderdale, W. J.; Gwaltney, S. R.; Beck, S.; Balkova, A.; Bernholdt, D. E.; Baeck, K.-K.; Rozyczko, P.; Sekino, H.; Huber, C.; Bartlett, R. J. ACES II is a program product of the Quantum Theory project, University of Florida. Integral packages included are VMOL (J. Almlof and P. R. Taylor), VPROPS (P.R. Taylor), ABACUS , (T. Helgaker, H. J. Aa Jensen, P. Jorgensen and P. R. Taylor).
Bofill, J. M.; Pulay, P. J. Chem. Phys.1989, 90, 3637–3646.
The statements, opinions and data contained in the journal International Journal of Molecular Sciences are solely
those of the individual authors and contributors and not of the publisher and the editor(s).
MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.
The statements, opinions and data contained in the journals are solely
those of the individual authors and contributors and not of the publisher and the editor(s).
MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.