Numerical and Theoretical Aspects of the DMRG-TCC Method Exemplified by the Nitrogen Dimer.

In this article, we investigate the numerical and theoretical aspects of the coupled-cluster method tailored by matrix-product states. We investigate formal properties of the used method, such as energy size consistency and the equivalence of linked and unlinked formulation. The existing mathematical analysis is here elaborated in a quantum chemical framework. In particular, we highlight the use of what we have defined as a complete active space-external space gap describing the basis splitting between the complete active space and the external part generalizing the concept of a HOMO-LUMO gap. Furthermore, the behavior of the energy error for an optimal basis splitting, i.e., an active space choice minimizing the density matrix renormalization group-tailored coupled-cluster singles doubles error, is discussed. We show numerical investigations on the robustness with respect to the bond dimensions of the single orbital entropy and the mutual information, which are quantities that are used to choose a complete active space. Moreover, the dependence of the ground-state energy error on the complete active space has been analyzed numerically in order to find an optimal split between the complete active space and external space by minimizing the density matrix renormalization group-tailored coupled-cluster error.

Introduction

The coupled-cluster (CC) theory has played a revolutionary role in establishing a new level of high accuracy in electronic structure calculations and quantum-chemical simulations. Despite the immense progress made in the field, computational schemes aiming at describing quasi-degenerate electronic structures of chemical systems are still unreliable. These multiconfiguration systems, also called strongly correlated systems, form one of the most challenging computational problems in quantum chemistry. Since these systems appear in various research areas, a reliable computational scheme is of major interest for the natural sciences. Recently, the computational advantages of the novel density matrix renormalization group-tailored coupled-cluster (DMRG-TCC) method restricted to single (S) and double (D) excitations were demonstrated on large and statically correlated systems by Veis et al.[1,2] Furthermore, computations have shown that the use of the DMRG-TCCSD method is indispensable for the DMRG in order to determine the proper structure of the low lying energy spectrum in strongly correlated systems.[2] In addition to these computational features, the DMRG-TCC approach is a promising candidate for a black-box quasi multireference scheme as considerable parts of the program already provide a routine procedure up to a numerical threshold. This would increase the accessibility for a broad class of researchers from various fields of study.

Although DMRG implementations already allow high precision multireference calculations on large complete active spaces (CAS), covering the majority of strongly correlated orbitals,[3] there is still a need for further analysis and developments in order to achieve a multireference routine procedure. In the setting of the DMRG-TCC method, a CAS DMRG solution is improved by means of an additional CC calculation performed on the remaining (external) orbital space. This CC correction, improving the description of dynamical correlation by the approximate solution, leaves the CAS part, i.e., the DMRG solution, invariant. We emphasize, however, that the presented implementation of the external CC corrections only take correlations with a single reference determinant into account, i.e., the considered CC amplitude equations are formulated with respect to one reference determinant. Despite the method’s dependence on the reference determinant (in its current version), we have noticed significant improvements for systems with multireference character (∼0.05Eh) via the CCSD correction on the external part compared to the single reference CCSD method on the full space. Nonetheless, the simplistic approach of the DMRG-TCC method to the multireference problem comes with a price. The DMRG-TCC, as a CAS method, does not correlate external amplitudes with the CAS amplitudes, i.e., contributions from the external part to excited determinants within the CAS are not present. Furthermore, in situations where the choice of a reference determinant becomes unclear, e.g., strong open-shell systems, the DMRG-TCC method could run into potential problems since it is based on a single reference formulation. Although the total spin can be fixed for the CAS part in the DMRG calculations (spin-adapted DMRG[4−7]), for the full orbital space it cannot be controlled through the external CC corrections presented in this work.

Common approaches to the strong correlation problem are provided by the multireference coupled-cluster (MRCC) theory based on the Jeziorski–Monkhorst ansatz.[8−10] The underlying idea of this ansatz is to include higher cluster excitations that are physically relevant but often more difficult to access in the usual single reference approach. To that end multiple determinants are employed in the reference state.[11] These multireference approaches can be roughly divided into three categories:[12] first, valence-universal approaches[13−21] (often also called genuine MRCC approaches); second, state-universal approaches;[22−27] and third, state-specific approaches.[28−49] Methods within the first two categories commonly suffer from so-called intruder states,[50−53] which leads to divergent behavior. Such methods furthermore require solving for a manifold of eigenstates, including several solutions that are irrelevant to the problem. These downsides can be overcome by state-specific approaches, however, they rely on an explicit inclusion of higher excitations. For a more detailed description of these active fields of research, we refer the reader to ref (9) and the references therein. An alternative multireference CC method that makes use of matrix product states and a modified DMRG algorithm is the linearized CCSD theory of Sharma and Alavi.[54] Furthermore, the pair CCD (pCCD) method—a CCD approach preserving electron pairs—besides being computationally inexpensive, can describe strong correlation, which the single reference CCD theory cannot. Nonetheless, the pCCD scheme lacks adequate dynamical correlation which was improved (by adding certain amplitudes) based on seniority of a determinant (number of unpaired electrons).[55] However, pairing merely the double excitations is not sufficient to describe the dissociation of the triple bond in the nitrogen dimer,[55] which is the content of this article. Higher order pairing schemes, however, allow a more effective treatment of strong correlations and are worth mentioning at this point.[56−62]

The mathematical analysis of CC schemes is far from being complete, especially with regard to multireference methods; however, many important steps have already been taken. The list of fundamental and mathematical chemistry articles aiming to describe the existence and nature of solutions of CC equations is too long to be summarized here. We will limit our discussion to a short selection of publications addressing such fundamental issues of chemistry.

As a system of polynomial equations the CC equations can have real or, if the cluster operator is truncated, complex solutions.[63,64] A standard tool to compute a solution of these nonlinear equations is the Newton–Raphson and the quasi-Newton method. However, these methods may diverge if the Jacobian or the approximated Jacobian become singular.[65] This is in particular the case when strongly correlated systems are considered. These and other related aspects of the CC theory have been addressed by Živković and Monkhorst[63,66] and Piecuch et al.[64] Significant advances in the understanding of the nature of multiple solutions of single-reference CC have been made by Kowalski and Jankowski[67] and by Piecuch and Kowalski.[68] An interesting attempt to address the existence of a cluster operator and cluster expansion in the open-shell case was done by Jeziorski and Paldus.[69]

The first advances within the rigorous realm of local functional analysis were performed by Schneider and Rohwedder, providing analyses of the closed-shell CC method for nonmulticonfiguration systems.[70−72] Since then, the local analysis of CC schemes was extended by Laestadius and Kvaal analyzing the extended CC method[73] and revisiting the bivariational approach[74,75] and Faulstich et al. providing the first local mathematical analysis of a multireference scheme, namely the CC method tailored by tensor network states (TNS-TCC).[76] As this mathematical branch of CC theory is very young, channeling abstract mathematical results into the quantum-chemistry community is a highly interdisciplinary task merging both fields. A first attempt in this direction was done by Laestadius and Faulstich linking the physical assumption of a HOMO–LUMO gap to the somewhat abstract Gårding inequality and in that context presenting key aspects of refs (70−73) from a more quantum chemical point of view.[77]

With this article, we aim to bridge the mathematical results in ref (76) of the TNS-TCC method into the quantum-chemistry community and extend these results with a numerical study on the complete active space dependence of the DMRG-TCCSD error. Furthermore, we derive formal properties of the TCC method.

DMRG-TCC Method

As a post-Hartree–Fock method, the TCC approach was introduced by Kinoshita et al. in ref (28) as an alternative to other multireference methods. It divides the cluster operator into a complete active space part, denoted Ŝ, and an external (ext) part T̂, i.e., the wave function is parametrized asSeparating the cluster operator into several parts goes back to Piecuch et al.[78,79] Note that the operators Ŝ, T̂ commute since this separation is merely a partition of the overall cluster operator. In this formulation the linked CC equations are given byComputing |ΨCAS⟩ = e|ΨHF⟩ first and keeping it fixed for the dynamical correction via the CCSD method restricts the above equations to |Ψμ⟩ not in the CAS, i.e., ⟨Ψ|Ψμ⟩ = 0 for all |Ψ⟩ in the CAS (we say that |Ψμ⟩ is in the L2-orthogonal complement of the CAS). We emphasize that this includes mixed states, e.g., |Ψ⟩ where |Ψ⟩ is an element of the CAS but |Ψ⟩ is not. We consider a CAS of N-electron Slater-determinants formed from the set of spin-orbitals . This is, in the mathematical sense, a subspace of the full configuration interaction (FCI) space, i.e., the space of all N-electron Slater-determinants formed from the entire set of spin-orbitals . We here assume the spin-orbitals to be eigenfunctions of the system’s Fock operator. Note that the following analysis can be applied to any single-particle operator fulfilling the properties used in ref (76)—not only the Fock operator. This mathematical analysis, among other things, rests on the structure of a one-particle operator with a distinct (and furthermore steerable) CAS-ext gap. As described below (in connection to Assumption A), choosing the Fock operator might lead to the inclusion of diffuse functions in the CAS.

Based on the single reference approach, the TCC method needs a large CAS to cover most of the static correlations. Since the size of the CAS scales exponentially with respect to the number of particles N, i.e., (for more details we refer the reader to ref (70)), an efficient approximation scheme for strongly correlated systems is indispensable for the TCC method to have practical significance. One of the most efficient schemes for static correlation is the DMRG method.[80] Going back to the physicists White and Martin,[3] it was introduced to quantum chemistry as an alternative to the CI or CC approach. However, the major disadvantage of the DMRG is that in order to compute dynamical correlation high bond dimensions (tensor ranks) may be necessary, making the DMRG a potentially costly method.[2,80] Nevertheless, as a TNS-TCC method, the DMRG-TCC approach is an efficient method since the CAS is supposed to cover the statically correlated orbitals of the system. This avoids the DMRG method’s weak point and allows to employ a much larger CAS compared to the traditional CI-TCC method. We remark here that some terminology has different meaning in mathematics, physics, and chemistry. The number of legs of a tensor is called the order of the tensor in mathematics, while it is called the rank of the tensor in physics. The rank of the matrix corresponds to the number of nonzero singular values after matricization in mathematics, i.e., the Schmidt number in physics.

A notable advantage of the TCC approach over some MRCC methods is that all excitation operators commute by construction. This is due to the fact that the Hartee–Fock method yields a single reference solution |ΨHF⟩, which implies that separating the cluster operator corresponds to a partition of excitation operators. Hence, Ŝ and T̂ commute. This makes the DMRG-TCC method’s analysis much more accessible than internally contracted MRCC methods and therewith facilitates establishing sound mathematical results.[76] We remark, however, that the computationally most demanding step of the DMRG-TCC calculation is the DMRG part, and its cost increases rapidly with k. Alternative to the dynamical correction via the CC approach, the DMRG-MRCI method in ref (81) utilizes an internally contracted CI algorithm different from a conventional CI calculation.

Formal Properties of the DMRG-TCC Method

It is desired that quantum-chemical computations possess certain features representing the system’s chemical and physical behavior. Despite their similarity, the CC and TCC method have essentially different properties, which are here elaborated. A basic property of the CC method is the equivalence of linked and unlinked CC equations. We point out that this equivalence is in general not true for the DMRG-TCCSD scheme. This is a consequence of the CAS ansatz since it yields mixed states, i.e., two particle excitations with one excitation into the CAS. The respective overlap integrals in the unlinked CC equations will then not vanish unless the single excitation amplitudes are equal to zero. Generalizing this result for rank complete truncations of order n we find that all excitation amplitudes need to be zero but for the nth one. This is somewhat surprising as the equivalence of linked and unlinked CC equations holds for rank complete truncations of the single-reference CC method.

For the sake of simplicity we show this results for the DMRG-TCCSD method. The general case can be proven in similar a fashion. We define the matrix representation T with elements Tμ,ν = ⟨Ψμ|e|Ψν⟩ for μ, ν ∉ CAS. Note that, as T̂ increases the excitation rank, T is an atomic lower triangular matrix and therefore not singular. Assuming that the linked CC equations hold, the nonsingularity of T yieldsAs the full projection manifold is complete under de-excitation, we obtain thatNote that the first term on the r.h.s. in eq together with the Hartree–Fock contribution from the sum, i.e., E0⟨Ψμ|e|ΨHF⟩, describe the unlinked CC equations. To analyze the remaining terms on the r.h.s. in eq we expand the inner products, i.e.,The first term in this expansion vanishes due to orthogonality. The same holds true for all terms where T̂ enters to the power of two or higher since an excitation of order two or higher acting on an at least singly excited Slater-determinant |Ψγ⟩ yields an at least 3-fold excited Slater-determinant. However, as the external space contains mixed states, we find that ⟨Ψμ|T̂|Ψγ⟩ is not necessarily zero, namely, for ⟨Ψμ| = ⟨Ψα| ∧ ⟨Ψβ| and |Ψγ⟩ = |Ψβ⟩ with α ∈ ext and β ∈ CAS. This proves the claim.

Subsequently, we elaborate the size consistency of the DMRG-TCCSD method. Let two DMRG-TCCSD wave functions for the individual subsystems A and B beThe corresponding energies are given byand the amplitudes fulfillin terms of the effective, similarity-transformed HamiltoniansThe Hamiltonian of the compound system of the noninteracting subsystems can be written as Ĥ = Ĥ + Ĥ. Since the TCC approach corresponds to a partitioning of the cluster amplitudes we note that forWith |ΨHF(⟩ = |ΨHF(⟩ ∧ |ΨHF(⟩, the energy of the compound systems can be written asIt remains to show thatsolves the Schrödinger equation, i.e., for all ⟨Ψμ(AB)|, it holds that . Splitting the argument into three cases, we note thatwhere ⟨Ψ(Ψ(| = ⟨Ψ(| ∧ ⟨Ψ(|. This proves the energy size consistency for the untruncated TCC method. From this we conclude the energy size consistency for the DMRG-TCCSD scheme, because the truncation only affects the product states ⟨Ψμ(Ψμ(| and these are zero in the above projection.

Looking at TCC energy expression we observe that due to the Slater–Condon rules, these equations are independent of CAS excitations higher than order three, i.e., amplitudes of Ŝ for n > 3. More precisely, due to the fact that in the TCCSD case external space amplitudes can at most contain one virtual orbital in the CAS, the TCCSD amplitude expressions become independent of Ŝ4, i.e.,where the primed variables a′, b′, c′, d′, e′ describe orbitals in the CAS, the nonprimed variable a describes an orbital in the external part and i, j, k, l, m, n are occupied orbitals. Note, this does not imply that we can restrict the CAS computation to a manifold characterizing excitations with rank less or equal to three as for strongly correlated systems these can still be relevant. However, it reduces the number of terms entering the DMRG-TCCSD energy computations significantly.

This work aligns with the originally introduced CI-TCCSD method taking only Ŝ for n = 1, 2 into account.[28] We emphasize that the additional consideration of Ŝ3 corresponds to an exact treatment of the CAS contributions to the energy. Furthermore, this consideration does not change the TCC method’s complexity, if the Ŝ3 amplitudes are available. This is due to the fact that including the CAS triple excitation amplitudes will not exceed the dominating complexities of the CCSD approach[82] nor of the DMRG method. However, the extraction of the CI-triples from the DMRG wave function is costly and a corresponding efficiency investigation is left for future work.

Analysis of the DMRG-TCC Method

In the sequel we discuss and elaborate mathematical properties of the TCC approach and their influence on the DMRG-TCC method. The presentation here is held brief and the interested reader is referred to ref (76) and the references therein for further mathematical details.

Complete Active Space Choice

As pointed out in the previous section, the TCC method relies on a well-chosen CAS, i.e., a large enough CAS that covers the system’s static correlation. Consequently, we require a quantitative measurement for the quality of the CAS, which presents the first obstacle for creating a nonempirical model since the chemical concept of correlation is not well-defined.[83] In the DMRG-TCC method, we use a quantum information theory approach to classify the spin-orbital correlation. This classification is based on the mutual-informationThis two particle entropy is defined via the von Neumann entropy S(ρ) = −Tr(ρ ln ρ) of the reduced density operators ρ{.[84] Note that the mutual-information describes two-particle correlations. For a more general connection between multiparticle correlations and ξ-correlations, we refer the reader to the work of Szalay et al.[84] We emphasize that in practice this is a basis dependent quantity, which is in agreement with the chemical definition of correlation concepts.[83] We identify pairs of spin-orbitals contributing to a high mutual information value as strongly correlated, the pairs contributing to the plateau region, i.e., a region in which the mutual information profile is constant, as nondynamically correlated and the pairs contributing to the mutual information tail as dynamically correlated (see Figure ). The mutual-information profile can be well approximated from a prior DMRG computation on the full system. Due to the size of the full system we only compute a DMRG solution of low bond dimension (also called tensor rank). These low-accuracy calculations, however, already provide a good qualitative entropy profile, i.e., the shapes of profiles obtained for low bond dimension, M, agree well with the ones obtained in the FCI limit. Here, we refer to Figures and 3 showing the single orbital entropy and mutual information profiles, respectively, for various M values and for three different geometries of the N2 molecule. The orbitals with large entropies can be identified from the low-M calculations providing a routine procedure to form the CAS including the strongly correlated orbitals.[85−87] In practice this is achieved by using the following dimension reduction protocol: We start with a very low bond dimension calculation carried out on the full orbital space. Based on the corresponding entropy profile and an a priori defined numerical threshold, a smaller set of orbitals is selected. In a subsequent step the same procedure is repeated on the reduced orbital set but with a larger bond dimension. This iterative dimension reduction protocol is a typical renormalization group based approach to refine the entropy spectrum that is also used in condensed matter physics.

Figure 3

Mutual information for r = 2.118a0 (blue), 2.700a0 (green), 3.600a0 (red) obtained for the full orbital space (k = 28) with DMRG for fixed M = 64, 256, 512 and for δεTr = 10–8, Mmax = 10 000.

Figure 2

Single orbital entropy for r = 2.118a0 (blue), 2.700a0 (green), 3.600a0 (red) obtained for the full orbital space (k = 28) with DMRG for fixed M = 64, 256, 512 and for δεTr = 10–8, Mmax = 10 000.

A central observation is that, for (i.e., k = N), the DMRG-TCCSD becomes the CCSD method and, for (i.e., k = K), it is the DMRG method. We recall that the CCSD method can not resolve static correlation and the DMRG method needs high tensor ranks for dynamically correlated systems. This suggests that the error obtains a minimum for some k with N ≤ k ≤ K, i.e., there exists an optimal choice of k determining the basis splitting and therewith the choice of the CAS. Note that this feature becomes important for large systems since high bond dimensions become simply impossible to compute with available methods.

Local Analysis of the DMRG-TCC Method

The CC method can be formulated as nonlinear Galerkin scheme,[70] which is a well-established framework in numerical analysis to convert the continuous Schrödinger equation to a discrete problem. For the DMRG-TCC method a first local analysis was performed in ref (76). There, a quantitative error estimate with respect to the basis truncation was established. Faulstich et al. showed under certain assumptions (Assumption A and B in the sequel) that the DMRG-TCC method possesses a locally unique and quasi-optimal solution (cf. section 4.1 in ref (76)). In case of the DMRG-TCC method the latter means: On a fixed CAS, the CC method tailored by a DMRG solution provides a truncation hierarchy that converges to the best possible dynamical correction to the given CAS. For a fixed basis set the CC solution tailored by a DMRG solution on a fixed CAS is up to a multiplicative constant the best possible solution in the approximation space defined by the basis set. In other words, the CC method provides the best possible dynamical correction for a given CAS solution such as a DMRG solution.

Note that local uniqueness ensures that for a fixed basis set, the computed DMRG-TCC solution is unique in a neighborhood around the exact solution. We emphasize that this result is derived under the assumption that the CAS solution is fixed. Consequently, for different CAS solutions we obtain in general different TCC solutions, i.e., different cluster amplitudes.

Subsequently, parts of the results in ref (76) are explained in a setting adapted to the theoretical chemistry perspective. The TCC function is given by f(t; s) = ⟨Ψμ|e–e–Ĥee|ΨHF⟩, for |Ψμ⟩ not in the CAS. Note that we use the convention where small letters s, t correspond to cluster amplitudes, whereas capital letters Ŝ, T̂ describe cluster operators. The corresponding TCC energy expression is given byConsequently, the linked TCC eqs then becomeWithin this framework the locally unique and quasi-optimal solutions of the TCC method were obtained under two assumptions (see Assumption A and B in ref (76)).

First, Assumption A requires that the Fock operator F̂ is bounded and satisfies a so-called . Note that spectral gap assumptions (cf. HOMO–LUMO gap) are standard in the analysis of dynamically correlated systems, and for a more detailed description of these properties in this context, we refer readers to ref (77). Second, in the theoretical framework[76] it is assumed that there exists a CAS-ext gap in the spectrum of the Fock operator, i.e., there is a gap between the kth and the k + 1st orbital energies. The CAS-ext gap (although in practice possibly very small) was sufficient for the analysis since the main purpose was to remove the HOMO–LUMO gap assumption and allow for quasi-degeneracy, which makes the general TCC approach applicable to multiconfiguration systems. Intuitively, this gap assumption means that the CAS captures the static correlation of the system.

However, in practice, an arbitrarily small gap is insufficient and needs to be complemented by a more detailed discussion (see Remark 10 in ref (76)). The crucial stability constant is not directly related to the CAS-ext gap ε – ε, nor to the HOMO–LUMO gap ε – ε. Due to the frozen CAS-amplitudes this stability constant becomes much larger and is roughly estimated by ε – ε. This improved stability provides accurate CC amplitudes and the improved gap is not destroyed e.g. by the existence of many diffuse functions around the LUMO-level (Fermi level). In this case, the CAS includes the diffuse functions. This might not be optimal but is the simplest choice and most importantly fulfills the stability condition. The issue of basis set optimization is discussed briefly in the conclusion but a more detailed discussion is left for future work.

Assumption B is concerned with the fluctuation operator Ŵ = Ĥ – F̂. This operator describes the difference of the Hamiltonian and a single particle operator, here chosen to be the Fock operator. Using the similarity transformed Ŵ and fixing the CAS amplitudes s, the mapis assumed to have a small enough Lipschitz-continuity constant (see eq 20 in ref (76)). The physical interpretation of this Lipschitz condition is at the moment unclear.

Error Estimates for the DMRG-TCC Method

A major difference between the CI and CC method is that the CC formalism is not variational. Hence, it is not evident that the CC energy error decays quadratically with respect to the error of the wave function or cluster amplitudes. Note that the TCC approach represents merely a partition of the cluster operator; however, its error analysis is more delicate than the traditional CC method’s analysis. The TCC-energy error is measured as a difference to the FCI energy. Let |Ψ*⟩ describe the FCI solution on the whole space, i.e., Ĥ |Ψ*⟩ = E|Ψ*⟩. Using the exponential parametrization and the above introduced separation of the cluster operator, we haveAn important observation is that the TCC approach ignores the coupling from the external space into the CAS. It follows that the FCI solution on the CAS |ΨFCICAS⟩ = exp(ŜFCI)|ΨHF⟩ is an approximation to the projection of |Ψ*⟩ onto the CASwhere P̂ = ∑μ∈CAS|Ψμ⟩⟨Ψμ| is the L2-orthogonal projection onto the CAS. For a reasonably sized CAS the FCI solution |ΨFCICAS⟩ is rarely computationally accessible and we introduce the DMRG solution on the CAS as an approximation of |ΨFCICAS⟩Tailoring the CC method with these different CAS solutions leads in general to different TCC solutions. In the case of |ΨFCICAS⟩, the TCC method yields the best possible solution with respect to the chosen CAS, i.e., f(tCC*;sFCI) = 0. This solution is in general different from tCC fulfilling f(tCC;sDMRG) = 0 and its truncated version tCCSD satisfying PGalf(tCCSD;sDMRG) = 0, where PGal denotes the l2-orthogonal projection onto the corresponding Galerkin space. In the context of the DMRG-TCC theory, the Galerkin space represents a truncation in the excitation rank of the cluster operator, e.g., DMRG-TCCD, DMRG-TCCSD, etc.

For the following argument, suppose that an appropriate CAS has been fixed. The total DMRG-TCC energy error ΔE can be estimated as[76]where each term on the r.h.s. in eq is now discussed. As a technical remark, the norms on either the Hilbert space of cluster amplitudes or wave functions are here simply denoted as ∥·∥. These norms are not just the l2- or L2-norm, respectively, but also measure the kinetic energy. It should be clear from context which Hilbert space is in question and we refer to ref (71) for formal definitions. The first term is defined aswhich describes the truncation error of the CCSD method tailored by |ΨDMRGCAS⟩. We emphasize that the dynamical corrections via the CCSD and the untruncated CC method are here tailored by the same CAS solution. Hence, the energy error Δε corresponds to a single reference CC energy error, which suggests an analysis similar to that of refs (70 and 72). Indeed, the Aubin–Nitsche duality method(88−90) yields a quadratic a priori error estimate in ∥tCCSD – tCC∥ (and in terms of the Lagrange mulitpliers; see Theorem 29 in ref (76)).

Second, we discuss the termHere, different CAS solutions with fixed external solutions are used to compute the energies. This suggests that ΔεCAS is connected with the errordescribing the approximation error of the DMRG solution on the CAS (see Lemma 27 in ref (76)). Indeedwith εμ = ε = ∑(λ – λ), for 1 ≤ n ≤ k, where λ are the orbital energies. The εμ are the (translated) Fock energies, more precisely, F̂ |Ψμ⟩ = (Λ0 + εμ)|Ψμ⟩, with Λ0 = ∑λ. Note that the wave function |ΨFCICAS⟩ is in general not an eigenfunction of Ĥ; however, it is an eigenfunction of the projected Hamiltonian P̂ĤP̂. Equation involves the exponential parametrization. This can be estimated by the energy error of the DMRG wave function, denoted , namelyIn section the energy error of the DMRG wave function is controlled by the threshold value δεTr, i.e., . Hence, for well chosen CAS the difference ∥|ΨDMRGCAS⟩ – |ΨFCICAS⟩∥ is sufficiently small such that holds. This again shows the importance of a well-chosen CAS. Furthermore, the last term in eq can be eliminated via orbital rotations, as it is a sum of single excitation amplitudes.

Finally, we considerSince (t*, s*) is a stationary point of we have . A calculation involving Taylor expanding around (t*, s*) (see Lemma 26 in ref (76)) yieldsNote that the above error is caused by the assumed basis splitting, namely, the correlation from the external part into the CAS is ignored. Therefore, the best possible solution for a given basis splitting (tCC*, sFCI) differs in general from the FCI solution (t*, s*).

Combining now the three quadratic bounds gives an overall quadratic a priori energy error estimate for the DMRG-TCC method. The interested reader is referred to ref (76) for a more detailed treatment of the above analysis.

On the k-Dependence of the Error Estimates

The error estimate outlined above is for a fixed CAS, i.e., a particular basis, splitting and bounds the energy error in terms of truncated amplitudes. Because the TCC solution depends strongly on the choice of the CAS, it is motivated to further investigate the k-dependence of the error ΔE. However, the above derived error bound has a highly complicated k-dependence since not only the amplitudes but also the implicit constants (in ≲) and norms depend on k. Therefore, the analysis in ref (76) is not directly applicable to take the full k-dependence into account.

In the limit where sDMRG → sFCI we obtain that tCC → tCC* since the TCC method is numerically stable, i.e., a small perturbation in s corresponds to a small perturbation in the solution t. Furthermore, if we assume that tCCSD ≈ tCC, which is reasonable for the equilibrium bond length of N2, the error can be bound asHere the subscript k on ΔE and C highlights the k-dependence. We remark that we here used the less accurate l2-structure on the amplitude space compared to the H1-structure in eq . This yields k-independent vectors (tCCSD, sDMRG) and (t*, s*), as well as an k-independent l2-norm. The k-depenence of C will be investigated numerically in more detail in section .

Splitting Error for N2

Including the k-dependence in the above performed error analysis explicitly is a highly nontrivial task involving many mathematical obstacles and is part of our current research. Therefore, we here extend the mathematical results from section with a numerical investigation on this k-dependence. Our study is presented for the N2 molecule using the cc-pVDZ basis, which is a common basis for benchmark computations developed by Dunning and co-workers.[91] Here we remark that in our calculations all electrons are correlated as opposed to the typical frozen-core calculation, where the two 1s orbitals are omitted from the full orbital space. We investigate three different geometries of the nitrogen dimer by stretching the molecule, thus the performance of DMRG-TCCSD method is assessed against DMRG and single reference CC methods for bond lengths r = 2.118a0, 2.700a0, and 3.600a0. In the equilibrium geometry the system is weakly correlated implying that single reference CC methods yield reliable results. For increasing bond length r the system shows multireference character, i.e., static correlations become more dominant. For r > 3.5a0 this results in the breakdown of single reference CC methods.[92] This breakdown can be overcome with the DMRG-TCCSD method once a large and well chosen CAS is formed, we therefore refer to the DMRG-TCCSD method as numerically stable with respect to the bond length along the potential energy surface (PES).

As mentioned before, the DMRG method is in general less efficient to recover dynamic correlations since it requires large computational resources. However, due to the specific CAS choice the computational resource for the DMRG part of the TCC scheme is expected to be significantly lower than a pure DMRG calculation for the same level of accuracy.

Computational Details

In practice, a routine application of the TCC method to strongly correlated molecular systems, i.e., to multireference problems, became possible only recently since it requires a very accurate solution in a large CAS including all static correlations. Tensor network state methods fulfill such a high accuracy criterion, but the efficiency of the TNS-TCCSD method strongly depends on various parameters of the involved algorithms. Some of these are defined rigorously while others are more heuristic from the mathematical point of view. In this section we present the optimization steps for the most important parameters of the DMRG-TCCSD method and outline how the numerical error study in section is performed.

As elaborated in sections and IV.A, the CAS choice is essential for the computational success of TNS-TCC methods. In addition, the error of the TNS method used to approximate the CAS part depends on various approximations. These include the proper choice of a finite dimensional basis to describe the chemical compound, the tensor network structure, and the mapping of the molecular orbitals onto the given network.[93] Fortunately, all these can be optimized by utilizing concepts of quantum information theory, introduced in section (see also the included references). In the following, we restrict the numerical study to the DMRG-TCCSD method but the results presented here should also hold for other TNS approaches.[93−97]

In the DMRG-TCCSD case, the tensor network topology in the CAS corresponds to a single branched tensor tree, i.e., a one-dimensional topology. Thus, permutations of orbitals along such an artificial chain effect the convergence for a given CAS choice.[98,99] This orbital-ordering optimization can be carried out based on spectral graph theory[100,101] by minimizing the entanglement distance,[102] defined as Idist = ∑I |i – j|2. In order to speed up the convergence of the DMRG procedure the configuration interaction based dynamically extended active space (CI-DEAS) method is applied.[93,99] In the course of these optimization steps, the single orbital entropy (S = S(ρ{)) and the two-orbital mutual information (I) are calculated iteratively until convergence is reached. The size of the active space is systematically increased by including orbitals with the largest single site entropy values, which at the same time correspond to orbitals contributing to the largest matrix elements of the mutual information. Thus, the decreasingly ordered values of S define the so-called CAS vector, which provides a guide in what order to extend the CAS by including additional orbitals. The bond dimensions M (tensor rank) in the DMRG method can be kept fixed or adapted dynamically (dynamic block state selection (DBSS) approach) in order to fulfill an a priori defined error margin.[103,104] Accurate extrapolation to the truncation free limit is possible as a function of the truncation error δεTr.[103,105]

In our DMRG implementation[106] we use a spatial orbital basis, i.e., the local tensor space of a single orbital is d = 4 dimensional. In this -representation an orbital can be empty, singly occupied with either a spin up or spin down electron, or doubly occupied with opposite spins. Note, in contrast to section we need N/2 spatial orbitals to describe an N-electron wave function and similar changes apply to the size of the basis set so that we use K ≡ K/2 from here on. The single orbital entropy therefore varies between 0 and ln d = ln 4, while the two-orbital mutual information varies between 0 and ln d2 = ln 16.

Next we provide a short description how to perform DMRG-TCCSD calculations in practice. Note that we leave the discussion on the optimal choice of k for the following sections.

First the CAS is formed from the full orbital space by setting k = K. DMRG calculations are performed iteratively with fixed low bond dimension (or with a large error margin) in order to determine the optimal ordering and the CAS vector as described above. Thus, the corresponding single-orbital entropy and mutual information are also calculated. These calculations already provide a good qualitative description of the entropy profiles with respect to the exact solution, i.e., strongly correlated orbitals can be identified.

Using a given N/2 < k < K we form the CAS from the Hartree–Fock orbitals and the first k – N/2 virtual orbitals from the CAS vector, i.e., orbitals with the largest single orbital entropy values. We emphasize that these orbitals contribute to the largest matrix elements in I. We carry out the orbital ordering optimization on the given CAS and perform a large-scale DMRG calculation with a low error threshold margin in order to get an accurate approximation of the |ΨFCICAS⟩. Note that the DMRG method yields a normalized wave function, i.e., the overlap with the reference determinant |ΨHF⟩ is not necessarily equal to one.

Using the matrix product state representation of |ΨDMRGFCI⟩ obtained by the DMRG method we determine the zero reference overlap, single and double CI coefficients of the full tensor representation of the wave function. Next, these are used to calculate the Ŝ1 and Ŝ2 amplitudes, which form the input of the forthcoming CCSD calculation.

In the following step the cluster amplitudes for the external part, i.e., T̂1 and T̂2, are calculated in the course of the DMRG-TCCSD scheme.

As we discus in the next section, finding the optimal CAS, i.e., k-splitting, is a highly nontrivial problem, and at the present stage we can only present a solution that is considered as a heuristic approach in terms of rigorous mathematics. In practice, we repeat steps 2–4 for a large DMRG-truncation error as a function of N/2 < k < K, thus we find local energy minima (see Figure ) using a relatively cheap DMRG-TCCSD scheme. Around such a local minimum we perform more accurate DMRG-TCCSD calculations by lowering the DMRG-truncation error in order to refine the optimal k. We also monitor the maximum number of DMRG block states required to reach the a priori defined DMRG-error margin as a function of k. Since it can happen that several k values lead to low error DMRG-TCCSD energies, while the computational effort increases significantly with increasing k we select the optimal k that leads to low DMRG-TCCSD energy but also minimizes the required DMRG block states. Using the optimal k value we perform large-scale DMRG-TCCSD calculation using a relatively tight error bound for the DMRG-truncation error.

Figure 4

Ground-state energy of the N2 molecule near the equilibrium geometry, r = 2.118a0, obtained with DMRG-TCCSD for 7 ≤ k ≤ 28 and for various DMRG truncation errors δεTr. The CCSD, CCSDT, and CCSDTQ reference energies are shown by dotted, dashed, and dashed–dotted lines, respectively. The CCSDTPQH energy (k = 28) is taken as a reference for the FCI energy. For δεTr = 10–5 the CAS was additionally formed by taking k orbitals according to increasing values of the single-orbital entropy values, i.e., inverse to the other CAS extensions. This is labeled by CAS↑ (see also Sec. .B.3).

We close this section with a brief summary of the numerically accessible error terms and relate them to equations presented in section . Note that the error analysis in section is presented for a given k, thus here the k dependence is also omitted.

For a given k split, the accuracy of |ΨDMRGCAS⟩ depends on the DMRG truncation error, δεTr. As has been shown in refs (103 and 105), the relative error, ΔErel = (EDMRG(CAS – EFCICAS)/EFCICAS is a linear function of δεTr on a logarithmic scale. Therefore, extrapolation to the FCI limit can be carried out as a function of δεTr. In addition, the error term appearing in eq can be controlled.

Note that terms appearing in eqs and 7 include FCI solutions of the considered system. However, for small enough and dynamically correlated systems, these FCI solutions can be well approximated. This is in particular the case for the nitrogen dimer near the equilibrium geometry with the here chosen basis set. The CI-coefficients are then extractable from the matrix product state representation of a wave function, e.g., |ΨDMRGCAS⟩ or |ΨFCICAS⟩. Note that calculating all CI-coefficients scales exponentially with the size of the CAS. However, since the system is dynamically correlated zeroth order, single and double excitation coefficients are sufficient. Hence, the error terms ∥|ΨDMRGCAS⟩ – |ΨFCICAS⟩∥ and ∥(ŜFCI – ŜDMRG()|ΨHF⟩∥ in eqs and 7, respectively, can be well approximated. We remark that this exponential scaling with the CAS size also effects the computational costs of the CAS CI-triples, which are needed for an exact treatment of the TCCSD energy equation. However, investigations of the influence of the CAS CI-triples on the computed energies are left for future work.

Results and Discussion

In this section, we investigate the overall error dependence of DMRG-TCCSD as a function of k and as a function of the DMRG-truncation error δεTr. For our numerical error study we perform steps 1–4 discussed in section for each N/2 < k < K. For each geometry r = 2.118a0, 2.700a0, and 3.600a0, we also carry out very high accuracy DMRG calculations on the full orbital space, i.e., by setting the truncation error to δεTr = 10–8 and k = K. This data is used as a reference for the FCI solution.

Entropy Study on the Full Orbital Space

We start our investigation by showing DMRG results for the full orbital space, i.e., the CAS is formed from k = K = 28 orbitals, and for various fixed M values and for δεTr = 10–8. In the latter case the maximum bond dimension was set to M = 10 000. In Figure a, we show the relative error of the ground-state energy as a function of the DMRG-truncation error on a logarithmic scale. For the FCI energy, EFCI, the CCSDTQPH reference energy is used given in ref (107). It is visible that the relative error is a linear function of the truncation error on a logarithmic scale, thus extrapolation to the truncation free solution can be carried out according to refs (103 and 105).

Figure 1

(a) Relative error of the ground-state energy as a function of the DMRG-truncation error on a logarithmic scale obtained for the full orbital space (k = K) with r = 2.118a0. (b) Maximum number of block states as a function of k for the a priori defined truncation error δεTr = 10–8 with r = 2.118a0 (blue), 2.700a0 (green), and 3.600a0 (red).

In Figures and 3, we present the sorted values of the single orbital entropy and of the mutual information obtained for fixed M = 64, 256, 512 and with δεTr = 10–8 for the three geometries. As can be seen in the figures, the entropy profiles obtained with low-rank DMRG calculations already resemble the main characteristics of the exact profile (M ≃ 10000). Therefore, orbitals with large single orbital entropies, also contributing to large matrix elements of I, can easily be identified from a low-rank computation. The ordered orbital indices define the CAS vector, and the CAS for the DMRG-TCCSD can be formed accordingly as discussed in section .

Taking a look at Figure , it becomes apparent that S shifts upward for increasing r indicating the higher contribution of static correlations for the stretched geometries. Similarly the first 50–100 matrix elements of I also take larger values for larger r while the exponential tail, corresponding to dynamic correlations, is less effected. The gap between large and small values of the orbital entropies gets larger and its position shifts rightward for larger r. Thus, for the stretched geometries more orbitals must be included in the CAS during the TCC scheme in order to determine the static correlations accurately. We remark here that the orbitals contributing to the high values of the single orbital entropy and mutual information matrix elements change for the different geometries according to chemical bond forming and breaking processes.[108]

Numerical Investigation of the Error’s k-Dependence

In order to obtain |Ψ*⟩ in the FCI limit, we perform high-accuracy DMRG calculations with δεTr = 10–8. The CAS was formed by including all Hartee–Fock orbitals and its size was increased systematically by including orbitals with the largest entropies according to the CAS vector. Orbitals with degenerate single orbital entropies, due to symmetry considerations, are added to the CAS at the same time. Thus, there are some missing k points in the following figures. For each restricted CAS we carry our the usual optimization steps of a DMRG scheme as discussed in section , with low bond dimension followed by a high-accuracy calculation with δεTr = 10–8 using eight sweeps.[93] Our DMRG ground-state energies for 7 < k < 28 together with the CCSD (corresponding to a DMRG-TCCSD calculation where k = N/2 = 7) and CCSDTQ reference energies, are shown in Figure near the equilibrium bond length, r = 2.118a0. The single-reference coupled cluster calculations were performed in NWChem,[109] we employed the cc-pVDZ basis set in the spherical representation. For k = K = 28 the CCSDTQPH energy was taken as a reference for the FCI energy.[107]

The DMRG energy starts from the Hartree–Fock energy for k = 7 and decreases monotonically with increasing k until the full orbital solution with k = 28 is reached. It is remarkable, however, that the DMRG-TCCSD energy is significantly below the CCSD energy for all CAS choices, even for a very small k = 9. The error, however, shows an irregular behavior taking small values for several different k values. This is due to the fact that the DMRG-TCCSD approach suffers from a methodological error, i.e., certain fraction of the correlations are lost, since the CAS is frozen in the CCSD correction. This supports the hypothesis of a k-dependent constant as discussed in section . Therefore, whether orbital k is part of the CAS or external part provides a different methodological error. This is clearly seen as the error increases between k = 10 and 15 although the CAS covers more of the system’s static correlation with increasing k. This is investigated in more detail in section .

Since several k-splits lead to small DMRG-TCCSD errors, the optimal k value from the computational point of view, is determined not only by the error minimum but also by the minimal computational time, i.e., we need to take the computational requirements of the DMRG into account. Note that the size of the DMRG block states contributes significantly to the computational cost of the DMRG calculation. The connection of the block size to the CAS choice is shown in Figure b, where the maximal number of DMRG block states is depicted as a function of k for the a priori defined truncation error margin δεTr = 10–8. Note that max(M) increases rapidly for 10 < k < 20. The optimal CAS is therefore chosen such that the DMRG block states are not too large and the DMRG-TCCSD provides a low error, i.e., is a local minimum in the residual with respect to k.

It is important to note that based on Figure the DMRG-TCCSD energy got very close to, or even dropped below, the CCSDT energy for several k values. Since close to the equilibrium geometry the wave function is dominated by a single reference character, it is expected that DMRG-TCCSD leads to even more robust improvements for the stretched geometries, i.e., when the multireference character of the wave function is more pronounced. Our results for the stretched geometries, r = 2.700a0 and 3.600a0, are shown in Figures , 3, 5, and 6. As mentioned in section , for larger r values static correlations gain importance signaled by the increase in the single orbital entropy in Figure . Thus, the multireference character of the wave function becomes apparent through the entropy profiles. According to Figure the DMRG-TCCSD energy for all k > 7 values is again below the CCSD computation and for k > 15 it is even below the CCSDT reference energy. For r = 3.600a0 the CC computation fluctuates with increasing excitation ranks and CCSDT is even far below the FCI reference energy, revealing the variational breakdown of the single-reference CC method for multireference problems. In contrast to this, the DMRG-TCCSD energy is again below the CCSD energy for all k > 7, but above the CCSDT energy. The error furthermore shows a local minimum around k = 19. For the stretched geometries static correlations are more pronounced, there are more orbitals with large entropies, thus the maximum number of DMRG block states increases more rapidly with k compared to the situation near the equilibrium geometry (see Figure b). Thus, obtaining an error margin within 1 μEh for k = 19 ≪ 28 leads to a significant save in computational time and resources. Here we remark that DMRG-TCCSD is a single-reference multireference method thus the choice of the reference determinant can effect its performance. In the our current study, however, we have verified that for d ≤ 4 and for all k values the weight of the Hartree–Fock determinant was significantly larger than all other determinants.

Figure 5

Ground-state energy of the N2 molecule with bond length r = 2.7a0, obtained with DMRG-TCCSD for 7 ≤ k ≤ 28 and for various DMRG truncation errors δεTr. The CCSD, CCSDT, and CCSDTQ reference energies are shown by dotted, dashed, and dashed–dotted lines, respectively. The DMRG energy with δεTr = 10–8 on the full space, i.e., k = 28, is taken as a reference for the FCI energy. For δεTr = 10–5, the CAS was additionally formed by taking k orbitals according to increasing values of the single-orbital entropy, i.e., inverse to the other CAS extensions. This is labeled by CAS↑ (see also section ).

Figure 6

Ground-state energy of the N2 molecule with bond length r = 3.6a0, obtained with DMRG-TCCSD for 7 ≤ k ≤ 28 and for various DMRG truncation errors δεTr. The CCSD, CCSDT, and CCSDTQ reference energies are shown by dotted, dashed, and dashed–dotted lines, respectively. The DMRG energy with δεTr = 10–8 on the full space, i.e., k = 28, is taken as a reference for the FCI energy. For δεTr = 10–5, the CAS was additionally formed by taking k orbitals according to increasing values of the single-orbital entropy, i.e., inverse to the other CAS extensions. This is labeled by CAS↑ (see also section ).

Effect of δεTr on the DMRG-TCCSD

In practice, we do not intend to carry out DMRG calculations in the FCI limit, thus usually a larger truncation error is used. Therefore, we have repeated our calculations for larger truncation errors in the range of 10–4 and 10–7. Our results are shown in Figures , 5, and 6. For small k the DMRG solution basically provides the Full-CI limit since the a priori set minimum number of block states Mmin ≃ 64 already leads to a very low truncation error. Therefore, the error of the DMRG-TCCSD is dominated by the methodological error. For k > 15 the effect of the DMRG truncation error becomes visible and for large k the overall error is basically determined by the DMRG solution. For larger δεTr between 10–4 and 10–5 the DMRG-TCCSD error shows a minimum with respect to k. This is exactly the expected trend, since the CCSD method fails to capture static correlation while DMRG requires large bond dimension to recover dynamic correlations, i.e., a low truncation error threshold. In addition, the error minima for different truncation error thresholds δεTr happen to be around the same k values. This has an important practical consequence: the optimal k-split can be determined by performing cheap DMRG-TCCSD calculations using large DMRG truncation error threshold as a function of k.

The figures furthermore indicate that ΔEGS has a high peak for 9 < k < 16. This can be explained by the splitting of the FCI space since this yields that the correlation from external orbitals with CAS orbitals is ignored. Thus, we also performed calculations for δεTr = 10–5 using a CAS formed by taking k orbitals according to increasing values of the single orbital entropy values in order to demonstrate the importance of the CAS extension. The corresponding error profile as a function of k near the equilibrium geometry is shown in Figure labeled by CAS↑. As expected, the improvement of DMRG-TCCSD is marginal compared to CCSD up to a very large k ≃ 23 split since ψDMRGCAS differs only marginally from ψHF.

Numerical Investigation on CAS-ext Correlations

Taking another look at Figure , we can confirm that already for small k values the most important orbitals, i.e., those with the largest entropies, are included in the CAS. In Figure , the sorted values of the mutual information obtained by DMRG(k) for 9 ≤ k ≤ 28 is shown on a semilogarithmic scale. It is apparent from the figure that the largest values of I change only slightly with increasing k, thus static correlations are basically included for all restricted CAS. The exponential tail of I corresponding to dynamic correlations, however, becomes more visible only for larger k values. We conclude, for a given k split the DMRG method computes the static correlations efficiently and the missing tail of the mutual information with respect to the full orbital space (k = 28) calculation is captured by the TCC scheme.

Figure 7

(a) Sorted values of the mutual information obtained by DMRG(k) for 9 ≤ k ≤ 28 on a semilogarithmic scale for N2 at r = 2.118a0. (b) Sorted 40 largest matrix elements of the mutual information obtained by DMRG(k) for 9 ≤ k ≤ 28 on a lin–lin scale for N2 at r = 2.118a0.

Correlations between the CAS and external parts can also be simulated by a DMRG calculation on the full orbital space using an orbital ordering according to the CAS vector. In this case, the DMRG left block can be considered as the CAS and the right block as the external part. For a pure target state, for example, the ground state, the correlations between the CAS and external part is measured by the block entropy, S(ρCAS() as a function of k. Here ρCAS( is formed by a partial trace on the external part of |ΨDMRGFCI⟩. The block entropy is shown in Figure a. The block entropy decays monotonically for k > 7, i.e, the correlations between the CAS and the external part vanish with increasing k. In contrast to this, when an ordering according to CAS↑ is used the correlation between CAS and external part remains always strong, i.e., some of the highly correlated orbitals are distributed among the CAS and the external part. Nevertheless, both curves are smooth and they cannot explain the error profile shown in Figure .

Figure 8

(a) Block entropy, S(ρCAS(), as a function of k for r = 2.118 ordering orbitals along the DMRG chain according to the same CAS and CAS↑ vectors as used in Figure . (b) e(k, δε) as a function of k of the nitrogen dimer near the equilibrium bond length for DMRG truncation error thresholds δε between 10–4 and 10–8.

Numerical Values for the Amplitude Error Analysis

Since correlation analysis based on the entropy functions cannot reveal the error profile shown in Figure , here we reinvestagate the error behavior as a function of N/2 ≤ k ≤ K but in terms of the CC amplitudes. Therefore, we also present a more detailed description of eq in section which includes the following terms:Here the valid index-pairs are μ = (, ), with = (i1, ..., i) ∈ {1, ..., N/2}, and = (a1, ..., a) ∈ {N/2 + 1, ..., K}. The excitation rank is given by |μ| = n where n = 1 stands for singles, n = 2 for doubles, and so on. The μ values are the labels of excitation operators τ̂ ≔ â†â, and τ̂ ≔ τ̂ ... τ̂. The corresponding amplitudes are given as t. For invalid index-pairs, i.e., index-pairs that are out of range, the amplitudes are always zero. The various amplitudes appearing in eq are calculated according to the following rules:

The tensor s*: amplitudes in the CAS(k) obtained by DMRG(δε* = 10–8) solution (represented by CI coefficients c*) for CAS(K)where i, i1, i2 ∈ {1, ..., N/2} and a, a1, a2 ∈ {N/2 + 1, ..., k}.

The tensor t*: amplitudes not in the CAS(k) obtained from the DMRG(δε* = 10–8) solution (represented by CI coefficients c*) for CAS(K) projected onto CAS(k), i.e., the complement (with respect to valid index-pairs) of s*.

The tensor sDMRG(k, δε): amplitudes in the CAS(k) are obtained by the DMRG(δε) solution (represented by CI coefficients c) for CAS(k). The amplitudes sDMRG(k, δε), sDMRG(k, δε) are the same as eq , but with c* → c, where i, i1, i2 ∈ {1, ..., N/2} and a, a1, a2 ∈ {N/2 + 1, ..., k}.

The tensor tCCSD(k, δε): amplitudes not in the CAS(k) obtained by TCCSD, i.e., the complement (with respect to valid index-pairs) of sDMRG(k, δε).

In Figure b we show the error e(k, δε) as a function of k of the nitrogen dimer near the equilibrium bond length. Note that the quantitative behavior is quite robust with respect to the bond dimension since the values only differ marginally. We emphasize that the error contribution in Figure is dominated by second term in eq since this is an order of magnitude larger than the contribution from the first and third terms in eq , respectively. The first term in eq is furthermore related to the usual T1 diagnostic in CC,[110] so it is not a surprise that a small value, ∼10–3, was found. Comparing this error profile to the one shown in Figure we can understand the irregular behavior and the peak in the error in ΔEGS between k = 9 and 17, and the other peaks for k > 17 but the error minimum found for k = 19 remains unexplained. Furthermore, we can conclude from Figure b that the quotient ΔEGS(k)/e(k, δε) is not constant. This indicates that the constants involved in section in particular the constant in eq in section is indeed k-dependent.

Conclusion

In this article we presented a fundamental study of the DMRG-TCCSD method. We showed that, unlike the single-reference CC method, the linked and unlinked DMRG-TCC equations are in general not equivalent. Furthermore, we showed energy size consistency of the TCC, DMRG-TCC, and DMRG-TCCSD method and gave a proof that CAS excitations higher than order three do not enter the TCC energy expression.

In addition to these computational properties of the DMRG-TCCSD method, we presented the mathematical error analysis performed in ref (76) from a quantum chemistry perspective. We showed local uniqueness and quasi optimality of DMRG-TCC solutions and highlighted the importance of the CAS-ext gap—a spectral gap assumption allowing to perform the analysis presented here. Furthermore, we presented a quadratic a priori error estimate for the DMRG-TCC method, which aligns the error behavior of the DMRG-TCC method with variational methods except for the upper bound condition. We emphasize that the DMRG-TCC solution depends strongly on the CAS choice. Throughout the analysis we neglected this dependence as we assumed an optimal CAS choice as indicated in section . The explicit consideration of this dependence in the performed error analysis carries many mathematical challenges, which are part of our current research. Therefore, we extended this work with a numerical study of the k-dependence of the DMRG-TCCSD error which showed also that the constants involved in the error estimation are most likely k-dependent. This stresses the importance of further mathematical work to include this dependence explicitly in the analysis.

We furthermore presented computational data of the single-site entropy and the mutual information that are used to choose the CAS. Our computations showed that these properties are qualitatively very robust, i.e., their qualitative behavior is well represented by means of a low-rank approximation, which is a computational benefit. The numerical investigation of the k-dependence of the DMRG-TCCSD error revealed that the predicted trend in section is correct. We demonstrated that the error indeed first decays (7 ≤ k ≤ 9) and then increases again (25 ≤ k ≤ 28) for low-rank approximations, i.e., 10–4 respectively 10–5. This aligns with the theoretical prediction based on the properties of the DMRG and single reference CC method. Additional to this general trend, the error shows oscillations. A first hypothesis is that this behavior is related to the ignored correlations in the transition k → k + 1. However, this was not able to be proven so far using entropy based measures but a similar irregular behavior can be detected by a cluster amplitude error analysis. Furthermore, such oscillations can be related to a bad reference function. Nonetheless, this scenario has here been ruled out since the Hartree–Fock determinant was found to be dominant in the CAS solution, i.e., the weight of the Hartree–Fock had largest weight in the CAS solution. The irregular behavior of the error minimum found for the DMRG-TCCSD method, therefore, could not be explained within this article and is left for future work. Despite the unknown reason for this behavior, we note that the error minima are fairly robust with respect to the bond dimension. Hence, the DMRG-TCCSD method can be extended with a screening process using low bond-dimension approximations to detect possible error minima.

On the other hand, an important feature that we would like to highlight here is that a small CAS (k = 9) yields a significant improvement of the energy and that the energies for all three geometries and all CAS choices outrun the single-reference CC method. In addition, the DMRG-TCCSD method avoids the breakdown of the CC approach even for multireference (strongly correlated) systems and, using concepts of quantum information theory, allows an efficient routine application of the method. Since the numerical error study showed a significant improvement for small CAS, we suspect the DMRG-TCCSD method to be of great use for larger systems with many strongly correlated orbitals as well as a many dynamically correlated orbitals.[1,2]

Finally, we remark, that besides the advantageous properties of the method there is a need for further analysis and developments in order to achieve our ultimate goal, i.e., to provide a black-box implementation of the DMRG-TCC method. Among these we highlight orbital rotations in the CAS through Fermionic mode transformation,[111] an automatic calculation of the best rank-1 representation of the DMRG wave function to be used as a reference state and the investigations of the influence of the CAS CI-triples on the computed energies. All these tasks are in progress.