Cluster perturbation theory. IV. Convergence of cluster perturbation series for energies and molecular properties

Pawłowski, Filip; Olsen, Jeppe; Jørgensen, Poul

doi:10.1063/1.5053622

The theoretical foundation has been developed for establishing whether cluster perturbation (CP) series for the energy, molecular properties, and excitation energies are convergent or divergent and for using a two-state model to describe the convergence rate and convergence patterns of the higher-order terms in the CP series. To establish whether the perturbation series are convergent or divergent, a fictitious system is introduced, for which the perturbation is multiplied by a complex scaling parameter z. The requirement for convergent perturbation series becomes that the energy or molecular property, including an excitation energy, for the fictitious system is an analytic, algebraic function of z that has no singularities when the norm |z| is smaller than one. Examples of CP series for the energy and molecular properties, including excitation energies, are also presented, and the two-state model is used for the interpretation of the convergence rate and the convergence patterns of the higher-order terms in these series. The calculations show that the perturbation series effectively become a two-state model at higher orders.

I. INTRODUCTION

In Paper I,¹ we introduced a new class of perturbation models—the cluster perturbation (CP) models—for which the major drawbacks of Møller-Plesset perturbation theory (MPPT)^2,3 and coupled cluster perturbation theory (CCPT)^4–6 have been overcome. The theoretical foundation for CP theory is given in Paper I.¹

In CP theory,¹ we consider a target excitation space relative to a Hartree-Fock (HF) state and partition the target excitation space into a parent and an auxiliary excitation space. The zeroth-order state in CP theory is a coupled cluster (CC) state in the parent excitation space, and we here assume that the target state is a CC state in the target excitation space. In CP theory, we determine perturbation series for the energy and for molecular properties, including excitation energies, in orders of the CC parent-state similarity-transformed fluctuation potential, where the zeroth-order term in the series is the energy or molecular property for the CC parent state and where the series formally converge to the energy or molecular property for the CC target state.

MPPT has recently been generalized to CCPT,^4–6 where the zeroth-order state is a CC state in the parent excitation space and where the target state is a CC state in the target excitation space. For both CCPT and CP theory, perturbation series for the energy are determined with the CC parent-state similarity-transformed fluctuation potential as the perturbation operator. However, in CCPT, terms are collected strictly as zeroth-order Fock operator contributions and first-order perturbation operator contributions, whereas in CP theory, a new, generalized order concept is introduced where one selected perturbation operator contribution is treated as a zeroth-order contribution. Using this generalized order concept, perturbation series can be determined on an equal footing for the energy and for molecular properties, contrary to CCPT for which perturbation series can only be determined for the ground-state energy.

In this paper, we examine the formal requirements for CP series to be convergent for the ground-state energy, for excitation energies, and for molecular properties. For the latter, we focus on first-order molecular properties and on molecular properties that can be described by the linear response function. We also discuss the theoretical foundation for setting up a two-state model^7,8 and for using this model to interpret the asymptotic convergence of CP perturbation series. The asymptotic convergence of CP series determines the convergence rate and the convergence patterns of the higher-order terms in the CP series. The two-state model has recently been thoroughly examined,⁹ and its convergence rate has been determined. Furthermore, the different archetypes that can arise for the higher-order terms in the perturbation series for different strengths of the interaction between the two states have been determined. We present numerical examples illustrating how the two-state model can be used to understand the convergence rate and the convergence pattern of the higher-order terms in CP series.

CP models¹ are characterized by the CC parent state, which is defined in the parent excitation space, and by an auxiliary excitation space. This can be expressed using a notation, where the parent excitation space is followed by the auxiliary space in parentheses. For example, CPSD(T) denotes a CP model with a coupled-cluster singles-and-doubles (CCSD) parent state and a triples auxiliary space. Furthermore, the notation CPSD(T) implies that a CC target state is used. If the auxiliary space is followed by a number, as for example in CPSD(T-3), the number denotes that perturbation corrections are determined through that order.

In Sec. II, we describe a general theoretical framework, which can be used for establishing whether the MPPT, CCPT, and CP energy series are convergent. The premises for setting up the two-state model for describing the asymptotic convergence of the perturbation series are also discussed. The requirements for energy series to be convergent, as well as the use of the two-state model for describing the asymptotic convergence of the energy series, have previously been considered for MPPT^7,8,10–12 and for CCPT.¹³

In this paper, we consider the requirements for CP series to be convergent not only for the energy, but also for excitation energies and molecular properties, where for the latter, we specifically consider first-order properties and properties that can be determined from the linear response function. To develop these requirements, we summarize in Sec. III how CP perturbation series is determined for the energy and for cluster amplitudes, and also for the Jacobian, by introducing the new, generalized order concept of CP.

In Sec. IV, the specific requirements for having convergent energy series in CP theory are developed. This is followed by numerical examples of CP series, for which the convergence rate and the convergence patterns of the higher-order terms are interpreted using the two-state model.⁹ In Sec. V, the requirements for convergent CP series are investigated for excitation energies and numerical examples are provided to illustrate how the asymptotic convergence of the series can be interpreted using the two-state model. Section VI contains a discussion of the requirements for convergent CP series for first-order molecular properties and for molecular properties that can be determined from the linear response function, along with numerical illustrations. Section VII contains a short summary and some concluding remarks.

II. PERTURBATION THEORY WITH A MØLLER-PLESSET PARTITIONING OF THE HAMILTONIAN

A. Requirements for convergent energy perturbation series

For a molecular system with the electronic Hamiltonian H, the electronic Schrödinger equation for the ground state |0⟩ can be written as

H | 0 (⟩ = E_{0} | 0 (⟩,

(1)

where E₀ is the ground-state energy. In MPPT, CCPT, and CP theory, the electronic Schrödinger equation is solved by projection using perturbation theory. The Møller-Plesset partitioning is used for the Hamiltonian,

H = f + Φ,

(2)

where f is the Fock operator and Φ is the fluctuation potential operator, and a perturbation series is determined for the ground-state energy,

E_{0} = E_{0}^{(0)} + \sum_{p = 1}^{\infty} E_{0}^{(p)},

(3)

where $E_{0}^{(0)}$ is the zeroth-order energy and $E_{0}^{(p)}$ is a term of order p in the perturbation. We discuss in this paper the theoretical foundation for establishing whether the perturbation series are convergent or divergent. Furthermore, we discuss the theoretical foundation for using the two-state model to describe the asymptotic convergence of the perturbation series.

We start by considering a fictitious system where a complex strength parameter z is multiplied on the fluctuation potential. The Hamiltonian for the fictitious system becomes

H (z) = f + z Φ .

(4)

The electronic Schrödinger equation for the fictitious system can be written as

H (z) | 0 (z) (⟩ = E_{0} (z) | 0 (z) (⟩,

(5)

where E₀(z) is an energy function in the complex plane z. For z = 1, the energy function gives the energy for the physical system,

E_{0} (1) = E_{0} .

(6)

For z = 0, we require that the energy function becomes the zeroth order energy for the physical system,

E_{0} (0) = E_{0}^{(0)} .

(7)

E₀(z) can be expanded in a Taylor series with z = 0 as the expansion point, giving

E_{0} (z) = E_{0}^{(0)} + \sum_{p = 1}^{\infty} z^{p} E_{0}^{(p)},

(8)

where $E_{0}^{(p)}$ is the pth-order expansion coefficient in the Taylor series. For z = 1, Eq. (8) gives the perturbation series in Eq. (3).

E₀(z) can have critical points z^c, at which E₀(z) is equal to another energy E_x(z),

E_{0} (z^{c}) = E_{x} (z^{c}) = E_{x 0} .

(9)

The critical point with the smallest norm |z^c| is called the primary critical point. The requirement for having a convergent ground-state energy series in Eq. (3) is that there are no critical points in the complex plane for E₀(z) inside the unit circle |z| ≤ 1. The primary critical point determines the convergence radius of the Taylor expansion in Eq. (8) and therefore has a central role for the determination of the asymptotic convergence of the perturbation series.

A state x with the energy E_x(z^c) that is degenerate with E₀(z^c) inside the unit circle |z| ≤ 1 is called an intruder state. If $Re (z^{c}) > 0$ ⁠, the state is called a front-door intruder state; conversely a back-door intruder state has $Re (z^{c}) < 0$ ⁠.⁸ Intruder states lead to divergent perturbation series. A state x with the energy E_x(z^c) that is equal to the ground-state energy, E₀(z^c), where |z^c| > 1, is called a crosser state. If a crosser state has $Re (z^{c}) \geq 0$ ⁠, it is called a front-door crosser state, while a crosser state with $Re (z^{c}) < 0$ is called a back-door crosser state.

In a finite dimensional space, the solution to the electronic Schrödinger equation for the fictitious system in Eq. (5) becomes a matrix equation,

H (z) C (z) = E (z) C (z),

(10)

where the Hamiltonian is partitioned into a zeroth-order Hamiltonian H₀ and a perturbation V,

H (z) = H_{0} + z V .

(11)

For configuration interaction (CI) target states as in MPPT the perturbation, V is a matrix representation of the fluctuation potential and it is symmetric, and E₀(z) is the lowest eigenvalue of the CI eigenvalue equation. The primary critical point is the point of degeneracy between the ground state and an excited state in the CI eigenvalue equation with the smallest norm |z^c|. For MPPT, the perturbation series in Eq. (3) thus describes a perturbation expansion of an eigenvalue equation.

For CCPT and CP theory, the target state is a CC state and E₀(z) is the CC energy, where amplitudes are determined by the set of non-linear cluster amplitude equations for the fictitious system. For CCPT, we have discussed in Ref. 13 the formal requirements for convergence of the CCPT energy series and have shown that intruder states occur as excitation operators, which give singularities within the unit circle for the perturbation-dependent Jacobian that is determined as the derivative of the CCPT cluster amplitude equations. In this paper, we discuss the theoretical foundation for having convergent ground-state CP series and show that the formal requirement for a convergent energy series is that the CP perturbation-dependent Jacobian J(z) does not have a singularity within the unit circle. We further discuss the requirements for having convergent CP series for molecular properties, including excitation energies. For CCPT and CP theory, the perturbation operator is the CC parent-state similarity-transformed fluctuation potential and the perturbation V in Eq. (11) is therefore non-symmetric.

B. The two-state model

The asymptotic convergence of a perturbation series determines the convergence rate of the perturbation series and the convergence patterns of the higher-order terms in the series. When studying the asymptotic convergence of MPPT series, it has been found that the correction vectors for these series for higher orders become nearly linearly dependent.^7,14 The asymptotic convergence of the MPPT series therefore can be studied using a two-state expansion. It has also been found that simple convergence patterns in MPPT and CCPT energy series can be described in terms of a two-state model^7,8 where the Hamiltonian matrix in Eq. (10) is set up for the fictitious system for the two states that compose the primary critical point. In this paper, we demonstrate that the asymptotic convergence of the CP energy series can be described in terms of a two-state model for the two states composing the primary critical point.

In this paper, we also consider the convergence of CP series for an excitation energy, ω_x = E_x − E₀, that is obtained by diagonalizing the Jacobian for the CC target state. We show that the requirement for a convergent CP excitation energy series for ω_x is that the excitation energy ω_x(z) = E_x(z) − E₀(z), which is obtained by diagonalizing the Jacobian for the fictitious system J(z), does not have for |z| < 1 another excitation energy, ω_y(z) = E_y(z) − E₀(z), that is degenerate with ω_x(z). The primary critical point for ω_x(z) is the point of degeneracy where ω_x(z^c) is equal to another excitation energy ω_y(z^c) with the smallest norm |z^c|. The primary critical point thus determines whether the CP series for the excitation energy ω_x is convergent.

For excitation energies, we also show in this paper that the asymptotic convergence of the CP series can be described in terms of a two-state model for the two states composing the primary critical point. The two-state model for excitation energies is obtained by setting up the perturbation-dependent Jacobian J(z) for the two states composing the primary critical point. We note that setting up the two-state model for the perturbation-dependent Jacobian is equivalent to setting up a two-state model for the equation-of-motion coupled cluster (EOM-CC) eigenvalue equation for the degenerate pair of excited state energies E_x(z) and E_y(z) since the Jacobian eigenvalue equation can be obtained from the EOM-CC eigenvalue equation¹⁵ by subtracting the ground-state energy E₀(z) from the diagonal.

To set up a two-state model that can be used for CP ground-state energy series and also for CP excitation energy series, we have to consider a more general two-state model than the one in Refs. 7 and 8, where the perturbation matrix in Eq. (11) was assumed to be symmetric. The Jacobian and the EOM-CC eigenvalue equations contain a non-symmetric matrix, and the interaction matrix between the two states composing the primary critical point therefore cannot be required to be symmetric. We have recently presented a general two-state model,⁹ where the interaction matrix is not symmetric and we summarize below the main features of this model.

1. The two-state model and its critical points

Using the same notation for the two-dimensional space, defined by the primary critical point, as for the full space in Eqs. (10) and (11), and assuming that the basis vectors of the two-dimensional space are orthonormal and diagonalize the zeroth-order Hamiltonian, the zeroth-order Hamiltonian and the perturbation can be written as⁹

H_{0} = (\begin{matrix} α & 0 \\ 0 & β + γ \end{matrix}),

(12)

V = (\begin{matrix} 0 & δ_{2} \\ δ_{1} & - γ \end{matrix}),

(13)

where α and β + γ are the two zeroth-order energies, and γ and δ₁, δ₂ are the gap shift and the coupling terms, respectively. The coupling terms, δ₁ and δ₂, can differ, and the perturbation matrix is therefore not required to be symmetric. We will in the following assume that β > α and β + γ > α, so the numerical order of the diagonal terms is the same in H₀ and H₀ + V.

The two eigenvalues of the matrix H₀ + zV are determined as^7,9

E_{\pm} (z) = \frac{α + β + (1 - z) γ}{2} \pm \frac{\sqrt{{(α - β - (1 - z) γ)}^{2} + 4 δ_{1} δ_{2} z^{2}}}{2} .

(14)

Equation (14) shows that the eigenvalues, as functions of z, depend on the product δ₁δ₂ and not on the individual coupling elements. We can therefore replace the individual coupling coefficients with a positive geometric average of these, δ, and a symmetry factor, σ,

δ = \sqrt{| δ_{1} δ_{2} |}

(15)

σ = \{\begin{matrix} + 1, & if δ_{1} δ_{2} \geq 0 \\ - 1, & if δ_{1} δ_{2} < 0, \end{matrix})

(16)

so

δ_{1} δ_{2} = σ δ^{2} .

(17)

The eigenvalues of Eq. (14) can be written in terms of σ and δ as

E_{\pm} (z) = \frac{α + β + (1 - z) γ}{2} \pm \frac{\sqrt{{(α - β - (1 - z) γ)}^{2} + 4 σ δ^{2} z^{2}}}{2} .

(18)

For the physical system, where z = 1, the square root in Eq. (18) becomes $\sqrt{{(α - β)}^{2} + 4 σ δ^{2}}$ and its deviation from |β − α| depends on the symmetry factor σ. Since the square root is a monotonically increasing function, it is first noted that $\sqrt{{(α - β)}^{2} + 4 σ δ^{2}}$ is larger than |β − α| for σ = 1 and smaller than |β − α| for σ = −1. For the symmetric perturbation, σ = 1, the lowest eigenvalue, E₋, is thus below α and the largest eigenvalue, E₊, is above β. For an asymmetric perturbation, σ = −1, and for $| δ | \leq \frac{| β - α |}{2}$ ⁠, one obtains a real lowest eigenvalue that is larger than α and a real largest eigenvalue that is lower than β, whereas for $| δ | > \frac{| β - α |}{2}$ ⁠, one obtains a pair of complex eigenvalues. For the ground-state energy, we have that the eigenvalue of interest is the lowest and that the perturbation expansion reduces the total energy. When applying the two-state model to the ground-state energy, we therefore only need to consider symmetric perturbations, σ = +1, and use Eq. (18) with this value of σ when deriving perturbation expansions for the two-state model. However, when applying the two-state model to CP series for excitation energies, we can encounter cases where the lower excitation energy increases and the higher excitation energy decreases when the perturbation is applied, and we then in addition also have to consider the asymmetric perturbations, σ = −1.

The primary critical points are defined by E₋(z) = E₊(z) and become

z_{\pm}^{c} = \frac{β - α + γ}{4 σ δ^{2} + γ^{2}} (γ \pm 2 \sqrt{- σ δ^{2}}) .

(19)

For symmetric perturbations, σ = 1, the critical point becomes a complex pair

z_{\pm}^{c} = \frac{β + γ - α}{4 δ^{2} + γ^{2}} (γ \pm 2 δ i),

(20)

whereas for asymmetric perturbations, σ = −1, the critical point becomes real

z_{\pm}^{c} = \frac{β - α + γ}{- 4 δ^{2} + γ^{2}} (γ \pm 2 δ) .

(21)

In this paper, we will report prototype examples of CP series for the ground-state energy and for excitation energies to illustrate that the asymptotic convergence of the series at high-orders effectively becomes a two-state problem. For the excitation energy series, we consider only examples where the asymptotic convergence of the CP series is described by symmetric perturbations. In the following, we therefore only consider the asymptotic convergence for symmetric perturbations and refer to Ref. 9 for a more complete treatment. The asymptotic rate of convergence, r, is equal to the inverse of the norm of the critical point, r = 1/|z^c|, so the location of the critical point does not only define whether the expansion is convergent or divergent, but it also defines the rate of convergence for a convergent expansion.

Since β + γ − α by assumption is positive, it is seen from Eq. (20) that the sign of the gap shift defines the position of the critical points in the complex plane: a positive gap shift leads to the critical point being located in the half plane with positive real values, whereas a negative gap shift leads to a critical point in the negative half plane.

2. Perturbation expansions of the energies for symmetric perturbations

In Ref. 9, we have identified perturbation expansions for the two-state model for both symmetric and non-symmetric perturbations. In this paper, we consider CP series only for symmetric perturbations, for which the energy corrections for E₋(1) can be written as^7,8

E_{-} (1) = \sum_{n = 0}^{\infty} E^{(n)},

(22)

E^{(0)} = α,

(23)

E^{(1)} = 0 .

(24)

The explicit expressions for E⁽ⁱ⁾, i = 2, 3, …, can be found in Ref. 8. To obtain a perturbation expansion for the excited state E₊(z), we see from Eq. (18) that the sum of the two eigenvalues contains only terms that are linear in the perturbation parameter,

E_{-} (z) + E_{+} (z) = α + β + (1 - z) γ .

(25)

Denoting the perturbation corrections of E₊(z) as $E_{+}^{(n)}$ and retaining the notation E⁽ⁿ⁾ for the perturbation expansion of E₋(z), we obtain the energy corrections for E₊(z) as⁹

E_{+}^{(0)} = β + γ,

(26)

E_{+}^{(1)} = - γ,

(27)

E_{+}^{(n)} = - E^{(n)}, n > 1 .

(28)

The second- and higher-order corrections for the excited eigenvalue are therefore equal to minus the correction E⁽ⁿ⁾ obtained for the lowest state.

3. Archetypes of convergence

In Ref. 9, we have analyzed the form of the energy corrections for the two-state perturbation expansion for general choices of γ, δ, and σ. For the symmetric perturbations, there are five archetypes of convergence patterns. Of these five archetypes—zigzag, interspersed zigzag, triadic, geometric, and ripples—all but the zigzag have been observed.⁹ We will examine the asymptotic convergence of CP series for the ground-state energy and for excitation energies, and we will as prototype examples consider the two archetypes: geometric and ripples. To have the background for examining the asymptotic convergence for these two archetypes, the properties and typical patterns of these two archetypes are given in Table I. For both archetypes, corrections and deviations are plotted for negative and positive gap shifts. The plots are logarithmic, and a simple color code is used to differentiate between positive (blue) and negative (red) corrections and deviations. To avoid unnecessary cluttering of the plots, the color codes are defined only on the first plot. The geometric and ripple archetypes arise when the absolute value of the gap shift γ is larger than the absolute value of the coupling element δ. If the gap shift is much larger than the coupling, the archetype is geometric, where the convergence exhibits a simple geometric form. For this archetype, a positive value of the gap shift implies that all corrections are negative, which in the table is denoted by (−), whereas a negative value of the gap shift leads to corrections of alternating sign, which is denoted by (1+, 1−). The sign of the gap shift can therefore be directly deduced from the signs of the corrections. In the ripple pattern, there are recurring ripples that are delineated by marked local minima in the size of the corrections and deviations. The number of orders spanned by a ripple, n^*, is proportional to the ratio $| \frac{γ}{δ} |$ ⁠. The signs of the corrections depend on the sign of the gap shift. If the gap shift is negative, the corrections have an alternating sign within a ripple, with the exception occurring at the boundary between two ripples, where two corrections have the same sign. If the gap shift is positive, all corrections have the same sign in a given ripple and the sign changes, when going from one ripple to the next.

TABLE I.

Archetypes of convergence patterns for two-state perturbation expansions (see text for details).

Archetype	Geometric	Ripples
Identification	\|γ\| ≫ \|δ\|	\|γ\| > \|δ\|
	Typical absolute corrections and deviations

Sign pattern	γ < 0: (1+, 1−)	γ < 0: (1+, 1−)
of corrections	γ > 0: (−)	γ > 0: (n^−, n^ +)
Period n^*	…	$2.5 + \sqrt{2} \|\frac{γ}{δ}\|$

Archetype	Geometric	Ripples
Identification	\|γ\| ≫ \|δ\|	\|γ\| > \|δ\|
	Typical absolute corrections and deviations

Sign pattern	γ < 0: (1+, 1−)	γ < 0: (1+, 1−)
of corrections	γ > 0: (−)	γ > 0: (n^−, n^ +)
Period n^*	…	$2.5 + \sqrt{2} \|\frac{γ}{δ}\|$

C. Outline

In the remainder of this paper, we will examine the formal requirements for having convergent CP series for the energy and also for molecular properties, including excitation energies. Furthermore, we will describe how the two-state model can be used to describe the asymptotic convergence of the CP series for the ground-state energy, excitation energies, and molecular properties. We also give numerical examples of CP series to demonstrate how the two-state model can be applied. However, before this can be done, we need to describe how perturbation series are obtained in CP theory for the energy and cluster amplitudes, and also for the Jacobian, and in particular how a new, generalized order concept is introduced in CP theory, compared to MPPT and CCPT, that allows perturbation series to be determined not only for the energy but also for molecular properties, including excitation energies.

III. CLUSTER PERTURBATION THEORY FOR THE ENERGY AND AMPLITUDE EQUATIONS AND FOR THE JACOBIAN

In CP theory, we consider a target excitation space, comprising all excitations through an excitation level t, that is partitioned into a parent excitation space, with excitations through a level p, and an auxiliary excitation space, with excitation levels p + 1 through t. The zeroth-order state in CP theory is a CC state in the parent excitation space. We assume here that the target state is a CC state^8,16 in the target excitation space.

The CC parent state can be written as

| {C C}^{*} (⟩ = e^{* T} | H F (⟩,

(29)

^{*} T =^{*} T_{1} + \dots +^{*} T_{p},

(30)

^{*} T_{i} = \sum_{μ_{i}} {^{*} t}_{μ_{i}} θ_{μ_{i}}, 1 \leq i \leq p .

(31)

The cluster operator ^*T_i contains the parent-state cluster amplitudes ${^{*} t}_{μ_{i}}$ and the many-body excitation operators $θ_{μ_{i}}$ that carry out excitations from the Hartree-Fock determinant |HF⟩ to excited determinants,

| μ_{i} (⟩ = θ_{μ_{i}} | H F (⟩ .

(32)

In Eqs. (31) and (32), i denotes an excitation level and μ_i an excitation at this level. The amplitudes of the CC parent state satisfy the cluster amplitude equations,

⟨) μ_{i} | e^{- * T} H_{0} e^{* T} | H F (⟩ = ⟨) μ_{i} | H_{0}^{^{*} T} | H F (⟩ = 0, 1 \leq i \leq p,

(33)

where the parent-state energy is

{}^{*}E_{0} = ⟨) H F | H_{0}^{^{*} T} | H F (⟩ .

(34)

In Eqs. (33) and (34), we have introduced the CC parent-state similarity-transformed Hamiltonian,

H_{0}^{^{*} T} = e^{-^{*} T} H_{0} e^{^{*} T} .

(35)

For p = 0, the parent excitation space is empty and ^*T vanishes, so the parent state becomes the Hartree-Fock state.

The CC target state in CP theory is parameterized using the CC parent state as the expansion point,

| C C (⟩ = e^{T} | H F (⟩ = e^{^{*} T + δ T} | H F (⟩ = e^{δ T} | {C C}^{*} (⟩,

(36)

where

T =^{*} T + δ T,

(37)

and

δ T = \sum_{i = 1}^{t} \sum_{μ_{i}} {δ t}_{μ_{i}} θ_{μ_{i}} .

(38)

The similarity-transformed Schrödinger equation for the CC target state can be written as

e^{- δ T} e^{-^{*} T} H_{0} e^{^{*} T} e^{δ T} | H F (⟩ = E_{0} | H F (⟩

(39)

and solved by projection in the target excitation space giving the cluster energy and amplitude equations,

E_{0} = ⟨) H F | e^{- δ T} e^{-^{*} T} H_{0} e^{^{*} T} e^{δ T} | H F (⟩ = ⟨) H F | e^{- δ T} H_{0}^{^{*} T} e^{δ T} | H F (⟩,

(40)

Ω_{μ_{i}} (δ t) = 0, 1 \leq i \leq t,

(41)

\begin{matrix} Ω_{μ_{i}} (δ t) & = & ⟨) μ_{i} | e^{- δ T} e^{-^{*} T} H_{0} e^{^{*} T} e^{δ T} | H F (⟩ \\ = & ⟨) μ_{i} | e^{- δ T} H_{0}^{^{*} T} e^{δ T} | H F (⟩, 1 \leq i \leq t . \end{matrix}

(42)

In CC theory, excitation energies are determined as eigenvalues of the CC response eigenvalue equation,¹⁷

J R_{x} = ω_{x} R_{x},

(43a)

L_{x} J = L_{x} ω_{x},

(43b)

L_{x} R_{y} = δ_{x y},

(43c)

where R_x and $L_{x}$ are right and left eigenvectors for an excited state x and ω_x is the excitation energy

ω_{x} = E_{x} - E_{0},

(44)

where E_x is the energy of excited state x. The left and right eigenvectors in Eq. (43) have been chosen to be biorthonormal. The Jacobian J in CP theory is parameterized with the CC parent state as an expansion point,

J_{μ_{i} ν_{j}} = ⟨) μ_{i} | [) e^{- δ T} H_{0}^{^{*} T} e^{δ T}, θ_{ν_{j}}]) | H F (⟩, i, j = 1,2, \dots, t,

(45)

where the cluster amplitudes of δT are determined from Eqs. (41) and (42). The Hamiltonian can be partitioned into a Fock operator, f, and a fluctuation potential operator, Φ, and $H_{0}^{^{*} T}$ in Eq. (35) can then be expressed as

H_{0}^{^{*} T} = f^{^{*} T} + Φ^{^{*} T} .

(46)

In CP theory, $Φ^{^{*} T}$ is used as the perturbation operator.¹

The extended parent-state Jacobian,

A_{μ_{i} ν_{j}} = ⟨) μ_{i} | [) H_{0}^{^{*} T}, θ_{ν_{j}}]) | H F (⟩, i, j = 1,2, \dots, t,

(47)

is a key quantity in CP theory. It appears in both a Baker-Campbell-Hausdorff (BCH) expansion of the cluster amplitude equations in Eqs. (41) and (42) and a BCH expansion of the Jacobian in Eq. (45). In CP theory, the extended parent-state Jacobian is partitioned into a zeroth-order and a first-order component,¹

A = J^{(0)} + J^{(1)},

(48)

where

\begin{matrix} J_{μ_{i} ν_{j}}^{(0)} & = & J_{μ_{i} ν_{j}}^{P} (1 - S_{i p}) (1 - S_{j p}) \\ + ⟨) μ_{i} | [) f^{^{*} T}, θ_{ν_{j}}]) | H F (⟩ δ_{μ_{i} ν_{j}} S_{i p} S_{j p}, 1 \leq i, j \leq t; p < t, \end{matrix}

(49)

\begin{matrix} J_{μ_{i} ν_{j}}^{(1)} & = & ⟨) μ_{i} | [) Φ^{^{*} T}, θ_{ν_{j}}]) | H F (⟩ (1 - S_{i p}) S_{j p} \\ + ⟨) μ_{i} | [) Φ^{^{*} T}, θ_{ν_{j}}]) | H F (⟩ S_{i p} (1 - S_{j p}) \\ + ⟨) μ_{i} | [) Φ^{^{*} T}, θ_{ν_{j}}]) | H F (⟩ S_{i p} S_{j p}, 1 \leq i, j \leq t; p < t, \end{matrix}

(50)

J_{μ_{i} ν_{j}}^{P} = ⟨) μ_{i} | [) H_{0}^{^{*} T}, θ_{ν_{j}}]) | H F (⟩, 1 \leq i, j \leq p,

(51)

and where we have introduced the integer step function S_ab,

S_{a b} = \{\begin{matrix} 0, & f o r a \leq b \\ 1, & f o r a > b . \end{matrix})

(52)

J^P is the CC parent state Jacobian, which constitutes the parent space component of J⁽⁰⁾ and contains a $Φ^{^{*} T}$ contribution.

In CP theory, the cluster amplitude equations in Eqs. (41) and (42) are solved order by order in $Φ^{^{*} T}$ with the exception that J^P, and thereby J⁽⁰⁾, is defined to be of zeroth order although it contains the parent-state similarity-transformed fluctuation potential projected onto the parent excitation space. It is the treatment of J^P as a zeroth-order contribution that is the key for the new order concept in CP theory, and it is this treatment of J^P that allows CP series to be determined not only for the ground-state energy but also for excitation energies and molecular properties.

To determine the CP series for the cluster amplitudes and for the Jacobian, we BCH expand the cluster amplitude equations in Eq. (42) and the Jacobian in Eq. (45) and use Eqs. (46)–(51) and then obtain

\begin{align} Ω_{μ_{i}} (δ t) & = \sum_{j = 1}^{t} \sum_{ν_{j}} J_{μ_{i} ν_{j}}^{(0)} {δ t}_{ν_{j}} + ⟨) μ_{i} | Φ^{^{*} T} | H F (⟩ S_{i p} \\ + \sum_{j = 1}^{t} \sum_{ν_{j}} J_{μ_{i} ν_{j}}^{(1)} {δ t}_{ν_{j}} + \frac{1}{2} ⟨) μ_{i} | [) [) Φ^{^{*} T}, δ T]), δ T]) | H F (⟩ \\ + \frac{1}{6} ⟨) μ_{i} | [) [) [) Φ^{^{*} T}, δ T]), δ T]), δ T]) | H F (⟩ \\ + \frac{1}{24} ⟨) μ_{i} | [) [) [) [) Φ^{^{*} T}, δ T]), δ T]), δ T]), δ T]) | H F (⟩, 1 \leq i \leq t, \end{align}

(53)

\begin{align} J_{μ_{i} ν_{j}} & = J_{μ_{i} ν_{j}}^{(0)} + J_{μ_{i} ν_{j}}^{(1)} + ⟨) μ_{i} | [) [) Φ^{^{*} T}, δ T]), θ_{ν_{j}}]) | H F (⟩ \\ + \frac{1}{2} ⟨) μ_{i} | [) [) [) Φ^{^{*} T}, δ T]), δ T]), θ_{ν_{j}}]) | H F (⟩ \\ + \frac{1}{6} ⟨) μ_{i} | [) [) [) [) Φ^{^{*} T}, δ T]), δ T]), δ T]), θ_{ν_{j}}]) | H F (⟩, 1 \leq i, j \leq t . \end{align}

(54)

A detailed derivation of Eqs. (53) and (54) can be found in Paper II.¹⁸

From the amplitude equations in Eqs. (41) and (53), we can determine a CP series in orders of $Φ^{^{*} T}$ for the cluster amplitudes,

{δ t}_{μ_{i}} = {δ t}_{μ_{i}}^{(1)} + {δ t}_{μ_{i}}^{(2)} + \dots,

(55)

with a vanishing zeroth-order term.¹ Substituting the cluster amplitude expansion in Eq. (55) into the energy in Eq. (40) and into the Jacobian in Eq. (54), we can determine CP series in orders of $Φ^{^{*} T}$ for the energy and the Jacobian,

E_{0} =^{*} E_{0} + \sum_{n = 1}^{\infty} E_{0}^{(n)},

(56)

J = J^{(0)} + J^{(1)} + J^{(2)} + \dots .

(57)

A detailed derivation of the CP series in Eqs. (55) and (56) can be found in Paper I¹ and of Eq. (57) in Paper II.¹⁸

In CP theory, the electronic Schrödinger equation for the CC target state in Eq. (39) is thus solved using projection giving the energy in Eq. (40) and the cluster amplitude equations in Eqs. (41) and (53). The cluster amplitude equations are solved using $Φ^{^{*} T}$ as the perturbation operator, where J⁽⁰⁾ in Eq. (49) is treated as a zeroth-order term. The energy in Eq. (40) can be written as

E_{0} (δ t) =^{*} E_{0} + Δ E_{0} (δ t),

(58)

where ^*E₀ is the parent state energy in Eq. (34) and ΔE₀(δt) is the energy correction that gives the CC target state energy,

Δ E_{0} (δ t) = ⟨) H F | Φ^{* T + δ T} - Φ^{* T} | H F (⟩ .

(59)

The parent state energy ^*E₀ is the zeroth-order energy and contains a fluctuation potential contribution. ΔE₀(δt) in Eq. (59) vanishes in the absence of the perturbation since ${δ t}_{μ_{i}}$ are then zero.

IV. CONVERGENCE OF THE CP ENERGY SERIES

We describe in Subsection IV A the theoretical foundation for examining whether the CP ground-state energy series are convergent or divergent and give in Subsection IV B numerical examples to illustrate the convergence of the CP ground-state energy series.

A. Theory

We now consider the fictitious system, described in Sec. II A, where we introduce a complex strength parameter z such that the perturbation becomes $z Φ^{^{*} T}$ ⁠. For the fictitious system, the electronic Schrödinger equation can be solved for the CC target state using projection, giving the energy and amplitude equations

E_{0} (δ t, z) =^{*} E_{0} (z) + Δ E_{0} (δ t, z),

(60)

Δ E_{0} (δ t, z) = z ⟨) H F | Φ^{* T + δ T} - Φ^{* T} | H F (⟩,

(61)

Ω_{μ_{i}} (δ t, z) = 0,

(62)

\begin{align} Ω_{μ_{i}} (δ t, z) & = \sum_{j = 1}^{t} \sum_{ν_{j}} J_{μ_{i} ν_{j}}^{(0)} (z) {δ t}_{ν_{j}} + ⟨) μ_{i} | z Φ^{^{*} T} | H F (⟩ S_{i p} \\ + z \sum_{j = 1}^{t} \sum_{ν_{j}} J_{μ_{i} ν_{j}}^{(1)} {δ t}_{ν_{j}} + \frac{1}{2} ⟨) μ_{i} | [) [) z Φ^{^{*} T}, δ T]), δ T]) | H F (⟩ \\ + \frac{1}{6} ⟨) μ_{i} | [) [) [) z Φ^{^{*} T}, δ T]), δ T]), δ T]) | H F (⟩ \\ + \frac{1}{24} ⟨) μ_{i} | [) [) [) [) z Φ^{^{*} T}, δ T]), δ T]), δ T]), δ T]) | H F (⟩, \end{align}

(63)

where the parent-state energy and amplitudes are determined from the equations

{}^{*}E_{0} (z) = ⟨) H F | e^{-^{*} T} H (z) e^{^{*} T} | H F (⟩,

(64)

⟨) μ_{i} | e^{-^{*} T} H (z) e^{^{*} T} | H F (⟩ = 0,

(65)

where H(z) is given in Eq. (4). J⁽⁰⁾ formally depends on the perturbation strength z through the z dependence of the perturbation operator contribution in the CC parent-state Jacobian. However, when the amplitude equations in Eqs. (62) and (63) are solved for ${δ t}_{μ_{i}}$ ⁠, this z dependence is treated explicitly and therefore does not need to be considered when the perturbation expansion of ${δ t}_{μ_{i}}$ is determined. For z = 1, the amplitude equations in Eqs. (62) and (63) become the amplitude equations for the physical system in Eqs. (41) and (53), and for z = 0, the amplitudes ${δ t}_{μ_{i}}$ vanish. For z = 1, the energy E₀(δt, z) becomes the energy in Eq. (58) of the physical system, and for z = 0, the energy E₀(δt, z) becomes ^*E₀.

E₀(δt, z) is an analytic algebraic function of the cluster amplitudes δt and the complex strength parameter z. The amplitudes δt depend on z, and we can therefore perform a Taylor series expansion of E₀(δt, z) with z = 0 as the expansion point,

E_{0} (δ t, z) =^{*} E_{0} + \sum_{n = 1}^{\infty} E_{0}^{(n)} z^{n} .

(66)

For z = 1, Eq. (66) becomes the CP series for the ground-state energy in Eq. (56). The Taylor series is convergent provided the derivative of the cluster amplitude equations with respect to the cluster amplitudes is non-singular inside the unit circle |z| ≤ 1.¹⁹ The derivative becomes

\begin{align} J_{μ_{i} ν_{j}} (z) & = \frac{d}{d t_{ν_{j}}} Ω_{μ_{i}} (δ t, z) = J_{μ_{i} ν_{j}}^{(0)} + z J_{μ_{i} ν_{j}}^{(1)} + ⟨) μ_{i} | [) [) z Φ^{^{*} T}, δ T]), θ_{ν_{j}}]) | H F (⟩ \\ + \frac{1}{2} ⟨) μ_{i} | [) [) [) z Φ^{^{*} T}, δ T]), δ T]), θ_{ν_{j}}]) | H F (⟩ \\ + \frac{1}{6} ⟨) μ_{i} | [) [) [) [) z Φ^{^{*} T}, δ T]), δ T]), δ T]), θ_{ν_{j}}]) | H F (⟩, \end{align}

(67)

where the z-dependence of J⁽⁰⁾ is suppressed since it is treated explicitly when the cluster amplitude equations are solved. For z = 1, Eq. (67) gives the CP Jacobian for the CC target state for the physical system in Eq. (54).

The complex Jacobian J(z) in Eq. (67) satisfies the eigenvalue equation

J (z) R_{x} (z) = ω_{x} (z) R_{x} (z),

(68a)

L_{x} (z) J (z) = L_{x} (z) ω_{x} (z),

(68b)

L_{x} (z) R_{x} (z) = 1,

(68c)

where

ω_{x} (z) = E_{x} (z) - E_{0} (z),

(69)

and E₀(z) is the energy of the CC target state in Eq. (60) and E_x(z) is the energy of an excited state x for the strength parameter z. The requirement for a convergent ground-state energy expansion is thus that the Jacobian J(z) has no singularities inside the unit circle |z| ≤ 1. At a singular point z^c of the Jacobian, the Jacobian has a vanishing eigenvalue,

ω_{x} (z^{c}) = E_{x} (z^{c}) - E_{0} (z^{c}) = 0,

(70)

and we have thus determined a critical point for E₀(δt, z). The determination of critical points for the energy function E₀(δt, z) is thus equivalent to determining singularities for the Jacobian J(z).

The search for the most important critical points for E₀(z) can be performed by searching for avoided crossings on the real axis $Re (z)$ for $E_{0} (Re (z))$ ⁠. When critical points are identified as singular points (zero eigenvalues) for the Jacobian J(z), the search for the avoided crossings for the primary critical point on the real axis $Re (z)$ for $E_{0} (Re (z))$ can be replaced by a search for the avoided crossings on the real axis $Re (z)$ for the smallest eigenvalue $ω_{x} (Re (z)) = E_{x} (Re (z)) - E_{0} (Re (z))$ of $J (Re (z))$ ⁠.

We now sketch how intruder states can arise when searching for avoided crossings on the real axis $Re (z)$ for the energy function $E_{0} (Re (z))$ ⁠. We will divide the search into two cases: a case where the parent space is either empty or contains singles, i.e., p = 0 or p = 1, and a case where p > 1.

For the first case with p = 0, Φ is the perturbation operator. The energy for the physical ground state, E₀(1), is the CP(SD…) energy, and the energy for the non-interacting system is E₀(0) = ∑_i $ε_{i}$ ⁠, where $ε_{i}$ denotes the orbital energy for an occupied orbital i. For p = 1, the singles cluster amplitudes in Eq. (33) vanish and Φ again becomes the perturbation operator. The energy for the physical ground state, E₀(1), is the CPS(DT…) energy, and the energy for the unperturbed system, E₀(0), is the Hartree-Fock energy. For both p = 0 and p = 1, the energy shift E₀(1) − E₀(0) is numerically large and negative. The energy shift E₀(1) − E₀(0) becomes, in particular, numerically large for electron-rich molecular systems. In Fig. 1(a), we have sketched a cartoon where we have marked the energy for the ground-state with a straight line that goes through E₀(0) and E₀(1), assuming linearity for the ground-state energy as a function of z. Furthermore, in Fig. 1(a), we have sketched a mechanism that gives a back-door intruder and a mechanism that gives a front-door intruder. For the mechanism that gives the back-door intruder, we consider a physical state E_y that is highly excited and diffuse. The energy of this state is nearly independent of the perturbation, and the energy of E_y therefore is nearly independent of z. Because of the large negative slope of the ground-state energy, the diffuse excited state with energy E_y gives rise to an avoided crossing at about z = −0.7 and thus to a back-door intruder. In Fig. 1(b), we have displayed E₀(δt, z) in the complex plane and marked with a circle the critical point associated with the back-door intruder. The circle has been located close to the real axis because of the very weak interaction between the highly excited diffuse state and the physical ground state. In Fig. 1(a), we have also sketched how a front-door intruder can arise as a result of an avoided crossing between a low lying physical state x with the energy E_x and the physical ground state with the energy E₀. Because of the strong interaction between these two states, the point of degeneracy representing the front-door intruder state in Fig. 1(b) is located away from the real axis but inside the unit circle resulting in a divergent CP series.

FIG. 1.

View large Download slide

For a parent excitation space, which contains at most single excitations, panel (a) sketches an energy diagram for a real perturbation strength $Re (z)$ that contains the ground-state energy E₀ and the excited-states energies E_x and E_y. The excited-state energy E_x has an avoided crossing with the ground-state energy E₀, which gives rise to a front-door intruder. For excited state E_y, the avoided crossing between E_y and E₀ gives rise to a back-door intruder. In panel (b), the cross denotes the critical point associated with the avoided crossing between E₀ and E_x, whereas the critical point associated with this back-door intruder is marked with a circle.

For p > 1, the cluster amplitudes of the CC parent state are non-vanishing and the perturbation operator becomes $Φ^{^{*} T}$ ⁠. The energy E₀(0) becomes the energy of the CC parent state, and E₀(1) becomes the energy of the CC target state. The energy difference E₀(1) − E₀(0) is small and negative, and in Fig. 2(a), we have sketched how a front-door intruder state can arise due to an avoided crossing between a low lying physical excited state E_x and the physical ground state E₀. Due to the strong interaction between the low lying excited state and the ground state, the point of degeneracy representing the front door intruder is located away from the real axis in Fig. 2(b) but inside the unit circle. When the p value of the parent space increases, the front-door intruder becomes better and better described by the CC parent state and the state giving rise to the front-door intruder eventually moves outside the unit circle and becomes a crosser state.

FIG. 2.

View large Download slide

For a parent excitation space, which contains at least single and double excitations, panel (a) sketches an energy diagram for real perturbation strength $Re (z)$ that contains the ground-state energy E₀ and the excited-state energy E_x. The excited-state energy E_x has an avoided crossing with the ground-state energy E₀ that gives rise to a front-door intruder. In panel (b), the cross denotes the critical point associated with the avoided crossing between E₀ and E_x.

For p > 1, the avoided crossing between the ground state and the highly diffuse excited state with energy E_y is displaced toward a large negative z value due to the numerically small negative slope of the ground-state energy E₀, and the back-door intruder state becomes a back-door crosser state. Back-door intruders therefore are much less frequent for parent excitation spaces that contain double and higher excitation levels. By increasing its diffuseness, the excited state will become less dependent of the perturbation and the back-door crosser state then will be displaced toward a numerically smaller negative z value, and the back-door intruder can eventually be reintroduced. The diffuseness of the excited state can be increased either by using a basis that contains more diffuse functions or by increasing the excitation level of the target excitation space.

Summarizing, the convergence of CP energy series is determined by states that are degenerate with the ground state in the complex plane for the energy in Eq. (60). Furthermore, these states can be determined by searching on the real axis $Re (z)$ for avoided crossings for the lowest eigenvalues of the Jacobian $J (Re (z))$ of Eq. (68). The asymptotic convergence of a CP energy series is determined by the point of degeneracy in the energy function in Eq. (60) with the smallest distance to the origin, i.e., the primary critical point. In Subsection IV B, we give numerical examples, which show how the asymptotic convergence of the CP energy series can be described using the two-state model for the fictitious system, which describes the interaction between the two states composing the avoided crossing associated with the primary critical point.

B. Numerical illustrations

We now present calculations of CP ground-state energy series to exemplify the above findings. We consider initially how back-door intruders can arise when the excitation level of the target excitation space is increased. We use the Ne atom in the aug-cc-pVDZ basis^20,21 as a first example and consider CP ground-state energy series where the parent excitation space contains only single excitations or both single and double excitations and where the target excitation space increases to contain up to hextuple excitations. The calculations were carried out with the LUCIA program²² where the generalized CC code²³ has been adapted to calculate CP ground-state and excitation energy series with general choices of the parent and the target excitation space. In Fig. 3, we have displayed logarithmic plots of the absolute errors of the ground-state energy for the series CPS(D-n), CPS(DT-n), …, CPS(DTQ56-n), where the errors are measured with respect to the CC target state energy. In Fig. 4, corresponding plots are given where the singles parent excitation space is replaced by a parent space containing both singles and doubles.

FIG. 3.

View large Download slide

Plots of the absolute errors for the ground-state energy for cluster perturbation series with a singles parent excitation space and an increasing size of the auxiliary excitation space for Ne in the aug-cc-pVDZ basis.

FIG. 4.

View large Download slide

Plots of the absolute errors for the ground-state energy for cluster perturbation series with a singles and doubles parent excitation space and an increasing size of the auxiliary excitation space for Ne in the aug-cc-pVDZ basis.

The plots in Figs. 3 and 4 show a geometric progression, which, according to Table I, is characteristic for a primary critical point, where |δ/γ| ≪ 1. The sign of the corrections alternates as exemplified by Fig. 5, where we have plotted the signed corrections for the convergent ground-state energy for the CPS(D-n) series and for the divergent CPS(DTQ56-n) ground-state energy series, which both exhibit an alternating pattern. Table I shows that the geometrical progressions in Figs. 3 and 4 can be associated with a negative gap shift and thereby with back-door intruder and crosser states.

FIG. 5.

View large Download slide

Signed energy corrections for cluster perturbation calculations for the ground-state energy series CPS(D-n) and CPS(DTQ56-n) for Ne in the aug-cc-pVDZ basis.

From Fig. 3, we see that the fastest convergence is obtained for the CPS(D-n) series. When higher excitation levels are considered in the target excitation space, the convergence becomes slower, and for the CPS(DTQ5-n) series, the ground state energy series diverges. The above trend is in accordance with the fact that the primary back-door crosser state in the energy function in Eq. (60) for the CPS(D-n) series shows up for a large negative z value and when the excitation level of the target excitation space increases, the diffuseness of the primary back-door crosser state increases and becomes more independent of the perturbation, and the negative z value for the avoided crossing becomes numerically smaller. For the CPS(DTQ5-n) series, the primary back-door crosser state becomes a back-door intruder state and the CP series diverges.

For the CPSD(T…) series in Fig. 4, the back-door intruder state shows up for hextuple excitations and thus for a higher excitation level than for the CPS(D…) series in Fig. 3. When increasing the excitation level of the parent excitation space to contain doubles, the intruder state can only have a very small doubles component and a higher target excitation level is therefore required to get a sufficiently diffuse excited state to introduce a back-door intruder. We also see from Figs. 3 and 4 that for a given target excitation space, the CP series converge faster for an increased excitation level of the parent excitation space. For example, the CPS(DT) series converges to a 10⁻⁶ a.u. accuracy at order 30, whereas for the CPSD(T) series, the same accuracy is obtained at order 8.

We now consider CP calculations of the ground-state energy series for BH at an internuclear distance of 3.5 a.u., i.e., about 1.5 times the equilibrium distance, using an aug-cc-pVDZ basis. The ground-state is multiconfigurational at this distance with a weight of 0.9 for the Hartree-Fock configuration 1σ²2σ²3σ² and a large 1σ²2σ²4σ² component. In Fig. 6, we have displayed logarithmic plots of the absolute errors for the ground-state energy for the CPS(DT-n) and CPSD(T-n) series. From Fig. 6, we see that both CP series converge fast and the energy function in Eq. (60) therefore has neither front-door nor back-door intruders. For FH molecule at 1.5 of its equilibrium distance and using the aug-cc-pVDZ basis, the CPS(DT-n) ground-state energy series diverges due to a back-door intruder. The multiconfigurational character is about the same for the FH and BH system. However, for the FH system, we have a large negative energy gap between the energy of the physical state, E₀(1), and the Hartree-Fock parent state energy, E₀(0), because FH is an electron-rich molecular system. This energy gap is numerically strongly reduced for the electron-poor BH system, and the intruder state for the FH system therefore becomes a crosser state for the BH system with $Re (z^{c}) < < - 1$ ⁠.

FIG. 6.

View large Download slide

Plots of the absolute errors for cluster perturbation calculations of the ground-state energy for the CPS(DT-n) and CPSD(T-n) series for BH at the internuclear distance 3.5 a.u. using the aug-cc-pVDZ basis.

The higher-order convergence for the CPS(DT-n) series in Fig. 6 is systematic, and from Table I, we recognize that it has a ripple pattern, for which |γ| > |δ|. For the plot of the CPS(DT-n) series in Fig. 6, we have introduced a color scheme for the absolute energy deviation, where a red circle indicates that the deviation is negative and a blue star denotes a positive energy deviation. We note that the deviations within a ripple have a constant sign. Comparing the deviation plot for the CPS(DT-n) series in Fig. 6 with the deviation plot for a ripple pattern with γ > 0 in Table I, the two plots show an amazing similarity. The ripple pattern is thus due to a primary front-door crosser state. The avoided crossing, which introduces this primary front-door crosser state, can be associated with the avoided crossing between the HF ground state and a singly excited state that has a large component of the 1σ²2σ²4σ² determinant. For the two states defining the avoided crossing, we know that |γ| > |δ|, where |δ| can be relatively large, since the interaction between these two states is large because BH has a multiconfigurational character at the distorted geometry.

The absolute error for the CPSD(T-n) ground-state energy is also given in Fig. 6 and shows a geometric progression. The color scheme for this curve shows that the energy deviations at each order have the same sign, and from Table I, we therefore recognize that the geometric progression is due to a primary front-door crosser state for which |γ| ≫ |δ|. For the CPSD(T-n) ground-state energy calculation, the avoided crossing for the primary critical point arises from the same two states that gave rise to the avoided crossing for the primary critical point for the CPS(DT-n) ground-state energy calculation. The only difference is that both the ground state and the excited state in the CPSD(T-n) calculation now are described to high accuracy in the parent excitation space and that the interaction between these two states therefore becomes very small. This means that |γ| ≫ |δ| and therefore leads to the geometric convergence pattern for the CPSD(T-n) series.

V. CONVERGENCE OF CP EXCITATION ENERGY SERIES

We describe in Subsection V A the theoretical foundation for examining the convergence of CP excitation energy series and give in Subsection V B numerical examples to illustrate the convergence of the CP series.

A. Theory

The Jacobian for the physical system in Eq. (57) is expanded in orders of the perturbation $Φ^{^{*} T}$ ⁠. The excitation energy ω_x and the right eigenvector R_x of the Jacobian eigenvalue equation in Eq. (43) can also be expanded in the orders of the perturbation,

ω_{x} = ω_{x}^{(0)} + ω_{x}^{(1)} + ω_{x}^{(2)} + \dots,

(71)

R_{x} = R_{x}^{(0)} + R_{x}^{(1)} + R_{x}^{(2)} + \dots,

(72)

and the CP perturbation series in Eqs. (71) and (72) can be determined by substituting Eqs. (53), (71), and (72) in Eq. (43) and solving Eq. (43) order by order in $Φ^{^{*} T}$ ⁠, assuming that R_x is intermediate normalized against $L_{x}^{(0)}$ ⁠. The zeroth-order equation reads

J^{(0)} R_{x}^{(0)} = ω_{x}^{(0)} R_{x}^{(0)},

(73a)

L_{x}^{(0)} J^{(0)} = L_{x}^{(0)} ω_{x}^{(0)},

(73b)

L_{x}^{(0)} R_{x}^{(0)} = 1,

(73c)

where J⁽⁰⁾ is given in Eq. (49). In a two-component form, referencing the parent (P) and auxiliary (A) excitation space components, we can write Eq. (73a) as

(\begin{matrix} J^{P} & 0 \\ 0 & ε_{A} \end{matrix}) (\begin{matrix} R_{x}^{P} \\ 0 \end{matrix}) = ω_{x}^{p} (\begin{matrix} R_{x}^{P} \\ 0 \end{matrix}),

(74)

where $ε_{A}$ is a diagonal matrix containing orbital energy differences. To obtain Eq. (74), we have used the CC parent state Jacobian eigenvalue equation,

J^{P} R_{x}^{P} = ω_{x}^{P} R_{x}^{P},

(75a)

L_{x}^{P} J^{P} = L_{x}^{P} ω_{x}^{P},

(75b)

L_{x}^{P} R_{x}^{P} = 1 .

(75c)

In zeroth order, we therefore have

ω_{x}^{(0)} = ω_{x}^{p},

(76)

R_{x}^{(0)} = (\begin{matrix} R_{x}^{P} \\ 0 \end{matrix}),

(77)

L_{x}^{(0)} = (\begin{matrix} L_{x}^{P} & 0 \end{matrix}) .

(78)

The CP series for the excitation energy ω_x in Eq. (71) is determined in Paper II.¹⁸

For the fictitious system, the Jacobian J(z) in Eq. (67) can be expanded in orders of the perturbation $z Φ^{^{*} T}$ and becomes

J (z) = J^{(0)} + z J^{(1)} + z^{2} J^{(2)} + \dots .

(79)

For z = 0, the Jacobian eigenvalue equation in Eq. (68) becomes the zeroth-order eigenvalue equation described by Eqs. (73)–(78), and for z = 1, the Jacobian eigenvalue equation in Eq. (68) becomes the Jacobian eigenvalue equation for the physical system in Eq. (43).

Let us assume that the Jacobian has no singularities inside the unit circle and that the CP expansion of δt in Eq. (55) therefore converges. ω_x(z) in Eq. (68) is then an analytic algebraic function in the complex plane of z and the eigenvector R_x(z). The eigenvalue ω_x(z) can be expanded in orders of z with z = 0 as the expansion point,

ω_{x} (z) = ω_{x}^{P} + \sum_{n = 1}^{\infty} w_{x}^{(n)} z^{n},

(80)

where $ω_{x}^{P}$ is an eigenvalue of J^P. For z = 1, Eq. (80) gives the CP expansion of the excitation energy in Eq. (71). The additional requirement for obtaining a convergent CP excitation energy expansion is therefore that there are no degeneracies of ω_x(z) inside the unit circle |z| ≤ 1 and that ω_x(z) therefore differs from the other excitation energies ω_y(z) inside the unit circle.

The search for degeneracies in the complex plane for ω_x(z) can in practice be performed searching for avoided crossings on the real axis $Re (z)$ for $J (Re (z))$ ⁠. Solving the Jacobian eigenvalue problem for $Re (z)$ is equivalent to determining the eigenvalues for the Hamiltonian matrix, but using an energy scale where the ground-state energy is subtracted, since the Jacobian eigenvalue equation can be obtained from the EOM-CC eigenvalue equation by subtracting the ground-state energy on the diagonal and removing the ground-state dimension.¹⁵ For z = 0, J(0) gives the parent state excitation energies, and for z = 1, J(1) gives the excitation energies of the physical system. The two-state model is for excitation energies obtained by carrying out a search for an avoided crossing between ω_x and another excitation energy ω_y and setting up the Jacobian eigenvalue equation for the fictitious system in the basis defined by the avoided crossing. We now consider the determination of the avoided crossing for ω_x and ω_y. Figure 7(a) depicts an energy diagram, where excitation energies for ω_x and ω_y are added to the ground-state reference energy E₀(z) = 0. Figure 7(a) sketches how an avoided crossing can lead to a front-door intruder between ω_x and ω_y. In Fig. 7(b), we have used a cross to mark the point of degeneracy in the complex plane that is associated with the avoided crossing for the excitation energy ω_x(R_x, z),

ω_{x} (z^{c}) = ω_{y} (z^{c}) = ω_{x y} .

(81)

For the avoided crossing, we have a large interaction between the two states and the primary critical point in Fig. 7(b) therefore has been given a larger imaginary component. We have also in Fig. 7(a) displayed the excitation energies for two other states, $ω_{v}$ and $ω_{w}$ ⁠. For these excitation energies, there is no point of degeneracies inside the unit circle. The point of degeneracy ω_xy for ω_x(z^c) and ω_y(z^c) thus does not affect whether the excitation energy series for $ω_{v}$ and $ω_{w}$ are convergent or divergent.

FIG. 7.

View large Download slide

Panel (a) sketches an energy diagram where excitation energies that are obtained by solving the Jacobian eigenvalue equation in Eq. (68) for real perturbation strength $Re (z)$ are plotted as a function of z. The excitation energies are plotted relatively to the ground-state energy E₀(z) = 0. Panel (a) shows an avoided crossing between ω_y and ω_x which introduces a front-door intruder for ω_x. In panel (b), the critical point associated with the avoided crossing is marked with a cross.

Figure 8(a) depicts an excitation energy diagram for excitation energy ω_x and ω_y, which leads to a back-door intruder. In contrast to Fig. 7(a), the degeneracy in the excitation energy is now due to a non-physical interaction as the strength parameter z has a negative sign.

FIG. 8.

View large Download slide

Panel (a) sketches an energy diagram where excitation energies that are obtained by solving the Jacobian eigenvalue equation in Eq. (68) for real perturbation strength $Re (z)$ are plotted as a function of z. The excitation energies are plotted relatively to the ground-state energy E₀(z) = 0. Panel (a) shows an avoided crossing between ω_y and ω_x which introduces a back-door intruder for ω_x. In panel (b), the critical point associated with the avoided crossing is marked with a cross.

In Figs. 7(a) and 8(a), the excitation ω_y introduces a front-door and a back-door intruder, respectively, in ω_x(R_x, z),

ω_{x} (z^{c}) = ω_{y} (z^{c}) = ω_{x y} .

(82)

The excitation energy ω_x will similarly give rise to a front-door and a back-door intruder in ω_y(R_y, z). The excitations composing the avoided crossing therefore are the same for ω_x(z) and for ω_y(z), and the two-state model therefore becomes identical for ω_x(z) and ω_y(z). This happens independently of whether the primary critical point is associated with an intruder or a crosser excitation. The target excitation energy for excitation x can be written as

ω_{x}^{T} = \sum_{j = 0}^{k} ω_{x}^{(j)} + \sum_{j = k + 1}^{\infty} ω_{x}^{(j)} = ω_{x k} + \sum_{j = k + 1}^{\infty} ω_{x}^{(j)},

(83)

where

ω_{x k} = \sum_{j = 0}^{k} ω_{x}^{(j)} .

(84)

Adding the target excitation energies for excitation x and y gives

ω_{x}^{T} + ω_{y}^{T} = ω_{x k} + ω_{y k} + \sum_{j = k + 1}^{\infty} () ω_{x}^{(j)} + ω_{y}^{(j)})) .

(85)

For the two-state model,⁹ we have for excitations x and y, in accordance with Eq. (28), that

{}^{T W O}ω_{x}^{(k)} = -^{T W O} ω_{y}^{(k)} .

(86)

Assuming that the two-state model can be applied to $ω_{x}^{T}$ and $ω_{y}^{T}$ starting from order k + 1, we have

{}^{T W O}ω_{x}^{(j)} = -^{T W O} ω_{y}^{(j)} = ω_{x}^{(j)} = - ω_{y}^{(j)}, j \geq k + 1,

(87)

and Eq. (85) then becomes

ω_{x}^{T} + ω_{y}^{T} = ω_{x k} + ω_{y k} .

(88)

Thus, even if the perturbation series for $ω_{x}^{T}$ and $ω_{y}^{T}$ are divergent, their sum for order k and larger will be equal to the sum of the two target excitation energies. In Subsection V B, we give numerical examples of excitation energies where we show that Eqs. (87) and (88) are satisfied to higher and higher accuracy when the order of the perturbation series increases, and this is also for divergent CP excitation energy series.

The CC parent and target state Jacobian J^P and J in Eqs. (75) and (43), respectively, and the perturbation scaled Jacobian J(z) of Eq. (68) are non-Hermitian, and we can therefore not guarantee that the Jacobian eigenvalue equations can be solved giving real eigenvalues and thus real excitation energies for truncated target excitation spaces. However, for the calculation of excitation energies for CC target states truncated at a given excitation level, the non-Hermiticity of the Jacobian have caused no problems with respect to obtaining real excitation energies and we will not discuss this issue in any further detail.

B. Numerical illustrations

To illustrate the above findings, we report CPS(D-n) calculations of the lowest excitation energies of ¹Π and ¹Σ⁺ symmetry for CO at the internuclear distance R_CO = 1.129 Å, using the aug-cc-pVDZ basis. The CPS(D-n) excitation energy corrections through an arbitrary order have been implemented using Psithon—a Python interface to the PSI4 program²⁴—and employing the Numerical Python Library (Numpy).²⁵ In Fig. 9, the absolute deviations are plotted for the two lowest excitations of ¹Π symmetry (A¹Π, E¹Π), and in Fig. 10, the absolute deviations are given for the four lowest excitations of ¹Σ⁺ symmetry (B¹Σ⁺, C¹Σ⁺, F¹Σ⁺, 4¹Σ⁺). From the figures, we see that for both the ¹Π and the ¹Σ⁺ excitations, the lowest-order corrections give a good estimate of the total excitation energy correction, which improves during the initial iterations.

FIG. 9.

View large Download slide

Absolute deviations from the target CCSD values for excitation energies (eV) of the two lowest-lying ¹Π states of the CO molecule, calculated using CPS(D-n) series through order 40 in aug-cc-pVDZ basis using a bond distance of 1.129 Å.

FIG. 10.

View large Download slide

Absolute deviations from the target CCSD values for excitation energies (eV) of the four lowest-lying ¹Σ⁺ states of the CO molecule, calculated using CPS(D-n) series through order 40 in aug-cc-pVDZ basis and a bond distance of 1.129 Å.

In Paper III,²⁶ we have benchmarked excitation energy corrections for the CPS(D-2) and CPS(D-3) models, recalling that the first-order excitation energy correction vanishes. For single-replacement dominated excitations the maximum and mean absolute errors (Δ_max, ${\bar{Δ}}_{a b s}$ ⁠), obtained by comparison with CCSD excitation energies, are (0.96 eV, 0.30 eV) for CPS(D-2) and (0.14 eV, 0.07 eV) for CPS(D-3). The errors for ¹Π and ¹Σ⁺ excitation energies are in line with the benchmark results. The CPS(D-3) model can therefore be used to get excitation energies of CCSD quality, in the sense that the difference between CPS(D-3) and CCSD excitation energies is of the same size or smaller than the effect of adding triples corrections to CCSD excitation energies.²⁶ Note that the CPS(D-2) model is identical to the configuration-interaction singles with a doubles correction [CIS(D)] model of Head-Gordon et al.²⁷

We now consider the higher-order convergence of the excitation energies and start with the excitations of ¹Π symmetry. After a monotonically convergent trend for initial orders, the A¹Π excitation shows a zig-zag pattern from orders 5-11 followed by a geometric convergence pattern. The excitation energy corrections have alternating sign, and according to Table I, the geometric progression therefore can be assigned to a back-door diffuse crosser state, where the distance of the critical point from the origin is large since the convergence for the geometric progression is fast. For the E¹Π excitation, the convergence shows a ripple pattern with a period of about 30. Furthermore, the signs of the corrections within a ripple alternate, and according to Table I, the convergence pattern for the E¹Π excitation can be assigned to a negative gap shift. Furthermore, as the period of the ripple is large, the size of the gap shift is much larger than the coupling.

Comparing the higher order convergence for the A ¹Π and E ¹Π excitations, we see from the curvature of the geometric convergence for the A ¹Π excitation that the convergence pattern actually is a ripple structure with a period that is larger than the period for the E ¹Π excitation. The ripple structure is due to an avoided crossing between the A ¹Π excitation and a diffuse highly excited state probably the same diffuse excited state that gave rise to the backdoor ripple structure for the E ¹Π excitation and with an avoided crossing located at a larger distance from the origin than the avoided crossing for the E ¹Π excitation.

We now consider the higher order convergence of the excitations of ¹Σ⁺ symmetry given in Fig. 10. After a converging trend for the first 3-4 orders, we see in Fig. 10 some orders with fluctuations in the excitation energy corrections. For the excitation pair (B¹Σ⁺, C¹Σ⁺), divergent ripple patterns of period 6 start at order 18, where the ripples are similar in the first period and become nearly identical in the later periods. For the excitation pair (F¹Σ⁺, 4¹Σ⁺), similar divergent ripple patterns start at order 11 with a period 10, and also for this pair, the patterns are similar and become nearly identical after the first period. According to the two-state model, ripple patterns are due to a significant coupling with intruder and crosser states, where the numerical value of the gap shift is larger than the interaction, i.e., |δ/γ| < 1 (see Table I).

To get a more detailed understanding of the divergence behavior for the two excitation pairs, we consider initially the excitation pair (B¹Σ⁺, C¹Σ⁺). At the coupled-cluster singles (CCS) level, the excitation energy to the B¹Σ⁺ state is 12.34 eV and the excitation vector has two dominant components: (5σ → 6σ) and (5σ → 7σ) with about equal weight. For the C¹Σ⁺ state, the CCS excitation energy is 12.80 eV and the excitation is described by the same two dominant excitations as for the B¹Σ⁺ state. At the CCSD level, the two excitations are separated, with the B¹Σ⁺ state having a CCSD excitation energy of 11.17 eV and the (5σ → 7σ) excitation as the dominant configuration, while the C¹Σ⁺ state has an excitation energy of 11.71 eV and the (5σ → 6σ) excitation as the dominant configuration. This leads to an avoided crossing with a front-door intruder for ω_B(z) and ω_C(z).

In Table II, convergence characteristics are given for the B¹Σ⁺ and the C¹Σ⁺ state through order 40. The two-state model predicts that the excitation energy corrections for the two states assigned to an avoided crossing become equal and of opposite sign. From Table II and Fig. 11, we see that the excitation energy corrections for the B¹Σ⁺ and the C¹Σ⁺ states for increasing orders become equal and have opposite sign. Furthermore, the two-state model also predicts that for a front-door intruder and the associated positive gap shift, the sign of the excitation energy corrections within a ripple is the same and that the sign changes when going from one ripple to the next. From Table II and Fig. 11, we see, as predicted, that the sign of the excitation energy corrections for higher orders changes with the period of 6 of the ripples.

TABLE II.

Excitation energy calculated through order k, ω_k = ω⁽⁰⁾ + ω⁽¹⁾ + ⋯ + ω^(k), for the two lowest ¹Σ⁺ states of CO molecule in aug-cc-pVDZ basis using a bond distance of 1.129 Å. Corrections at a given order, ω^(k), and the sum of excitation energies for both states at a given order, ω_k(B¹Σ⁺) + ω_k(C¹Σ⁺), are also reported. Note that ω⁽⁰⁾ = ω^CCS, ω⁽¹⁾ = 0, and ω₂ = ω^CIS(D). All results in eV.

	B¹Σ⁺		C¹Σ⁺
k	ω^(k)	ω_k	ω^(k)	ω_k	ω_k(B¹Σ⁺) + ω_k(C¹Σ⁺)
0	12.347	12.347	12.799	12.799	25.146
1	0.000	12.347	0.000	12.799	25.146
2	−0.898	11.448	−1.182	11.618	23.066
3	−0.131	11.317	−0.193	11.425	22.742
4	−0.110	11.208	0.189	11.614	22.821
5	0.032	11.239	0.064	11.677	22.917
6	−0.124	11.115	0.063	11.740	22.855
7	0.013	11.128	0.051	11.790	22.918
8	−0.056	11.072	−0.028	11.763	22.834
9	0.027	11.099	0.050	11.813	22.911
10	0.009	11.108	−0.086	11.727	22.834
11	0.051	11.159	0.022	11.749	22.908
12	0.020	11.179	−0.089	11.660	22.839
13	0.067	11.246	−0.004	11.656	22.902
14	−0.021	11.225	−0.037	11.619	22.844
15	0.053	11.279	−0.001	11.619	22.897
16	−0.073	11.205	0.026	11.644	22.849
17	0.007	11.212	0.036	11.680	22.892
18	−0.084	11.128	0.047	11.726	22.854
19	−0.037	11.092	0.069	11.796	22.887
20	−0.036	11.056	0.007	11.803	22.859
21	−0.027	11.029	0.051	11.854	22.883
22	0.037	11.067	−0.058	11.796	22.862
23	0.043	11.110	−0.026	11.770	22.880
24	0.072	11.182	−0.087	11.684	22.865
25	0.110	11.292	−0.098	11.585	22.877
26	0.032	11.324	−0.042	11.544	22.868
27	0.083	11.407	−0.075	11.469	22.875
28	−0.058	11.349	0.052	11.521	22.870
29	−0.056	11.293	0.059	11.580	22.873
30	−0.122	11.172	0.119	11.699	22.871
31	−0.189	10.983	0.190	11.889	22.872
32	−0.090	10.894	0.089	11.978	22.872
33	−0.147	10.747	0.146	12.125	22.871
34	0.041	10.788	−0.040	12.085	22.873
35	0.106	10.894	−0.108	11.977	22.871
36	0.185	11.078	−0.183	11.795	22.873
37	0.358	11.437	−0.361	11.434	22.871
38	0.207	11.643	−0.204	11.230	22.873
39	0.291	11.934	−0.293	10.937	22.871
40	0.021	11.955	−0.019	10.918	22.873
41	−0.182	11.773	0.179	11.097	22.871

	B¹Σ⁺		C¹Σ⁺
k	ω^(k)	ω_k	ω^(k)	ω_k	ω_k(B¹Σ⁺) + ω_k(C¹Σ⁺)
0	12.347	12.347	12.799	12.799	25.146
1	0.000	12.347	0.000	12.799	25.146
2	−0.898	11.448	−1.182	11.618	23.066
3	−0.131	11.317	−0.193	11.425	22.742
4	−0.110	11.208	0.189	11.614	22.821
5	0.032	11.239	0.064	11.677	22.917
6	−0.124	11.115	0.063	11.740	22.855
7	0.013	11.128	0.051	11.790	22.918
8	−0.056	11.072	−0.028	11.763	22.834
9	0.027	11.099	0.050	11.813	22.911
10	0.009	11.108	−0.086	11.727	22.834
11	0.051	11.159	0.022	11.749	22.908
12	0.020	11.179	−0.089	11.660	22.839
13	0.067	11.246	−0.004	11.656	22.902
14	−0.021	11.225	−0.037	11.619	22.844
15	0.053	11.279	−0.001	11.619	22.897
16	−0.073	11.205	0.026	11.644	22.849
17	0.007	11.212	0.036	11.680	22.892
18	−0.084	11.128	0.047	11.726	22.854
19	−0.037	11.092	0.069	11.796	22.887
20	−0.036	11.056	0.007	11.803	22.859
21	−0.027	11.029	0.051	11.854	22.883
22	0.037	11.067	−0.058	11.796	22.862
23	0.043	11.110	−0.026	11.770	22.880
24	0.072	11.182	−0.087	11.684	22.865
25	0.110	11.292	−0.098	11.585	22.877
26	0.032	11.324	−0.042	11.544	22.868
27	0.083	11.407	−0.075	11.469	22.875
28	−0.058	11.349	0.052	11.521	22.870
29	−0.056	11.293	0.059	11.580	22.873
30	−0.122	11.172	0.119	11.699	22.871
31	−0.189	10.983	0.190	11.889	22.872
32	−0.090	10.894	0.089	11.978	22.872
33	−0.147	10.747	0.146	12.125	22.871
34	0.041	10.788	−0.040	12.085	22.873
35	0.106	10.894	−0.108	11.977	22.871
36	0.185	11.078	−0.183	11.795	22.873
37	0.358	11.437	−0.361	11.434	22.871
38	0.207	11.643	−0.204	11.230	22.873
39	0.291	11.934	−0.293	10.937	22.871
40	0.021	11.955	−0.019	10.918	22.873
41	−0.182	11.773	0.179	11.097	22.871

FIG. 11.

View large Download slide

Excitation energy corrections, ω^(k) (eV), for the B¹Σ⁺ and C¹Σ⁺ states of the CO molecule calculated using CPS(D-n) series through order 40 in aug-cc-pVDZ basis using a bond distance of 1.129 Å. The horizontal line, where the excitation energy is zero, is a mirror plane for the higher-order excitation energy corrections.

Assuming that CP excitation energy series can be described asymptotically by the two-state model, the sum of the excitation energies ω_k(B¹Σ⁺) + ω_k(C¹Σ⁺) in the last column of Table II will converge toward the sum of the CCSD excitation energies,

ω_{B}^{C C S D} + ω_{C}^{C C S D} = 22.88 eV .

This sum is clearly approached for increasing k in the last column of Table II. The deviation of ω_k(B¹Σ⁺) + ω_k(C¹Σ⁺) from 22.88 eV for increasing k gives an indication when the asymptotic region is reached and shows that this happens at about order 20, where the ripples pattern in Fig. 10 for the two states becomes identical and where the excitation energy corrections for the B¹Σ⁺ and the C¹Σ⁺ state in Fig. 11 become mirror images of each other.

Although the divergence pattern for the (B¹Σ⁺, C¹Σ⁺) and the (F¹Σ⁺, 4¹Σ⁺) pairs looks similar, their origin differs. The divergence for the (B¹Σ⁺, C¹Σ⁺) pair is due to an avoided crossing between the B¹Σ⁺ and C¹Σ⁺ excitations that introduces a front-door intruder state, whereas the divergence for the (F¹Σ⁺, 4¹Σ⁺) pair is due to a back-door intruder. To confirm that this is the case, we report in Table III convergence characteristics for the (F¹Σ⁺, 4¹Σ⁺) excitation pair. The two-state model for a back door intruder predicts that the excitation energy corrections alternate within a ripple and that two corrections with the same sign occur at the border between two consecutive ripples. From Table III, we also see that for orders 12 and 13, the excitation energy corrections for both the F¹Σ⁺ and the 4¹Σ⁺ state have the same sign and are followed by 7 orders with an alternating sign, 2 orders with the same sign, 8 orders with an alternating sign, 2 orders with the same sign, and so on. This is precisely what the two-state model predicts for a back-door intruder state with a small interaction between the states. In Fig. 12, we have plotted excitation energy corrections for the F¹Σ⁺ and the 4¹Σ⁺ excitations. For orders larger than 10, the excitation energy corrections become mirror images and display the alternating trend described above. Comparing the higher-order convergence for the (B¹Σ⁺, C¹Σ⁺) and the (F¹Σ⁺, 4¹Σ⁺) pairs, we further observe that while the asymptotic convergence for the (B¹Σ⁺, C¹Σ⁺) pair starts around order 20, for the (F¹Σ⁺, 4¹Σ⁺) pair, it starts at around order 10. We also note that the ripples for the ripple pattern for the (F¹Σ⁺, 4¹Σ⁺) pair is right-tilted as this pattern is connected to a fast convergent perturbation series, whereas the ripple pattern for the B¹Σ⁺ and C¹Σ⁺ states is a protracted pattern with ripples that are nearly symmetric since the ripple patterns in these cases are connected to a slowly divergent perturbation series.

TABLE III.

Excitation energy calculated through order k, ω_k = ω⁽⁰⁾ + ω⁽¹⁾ + ⋯ + ω^(k), for the F¹Σ⁺ and 4¹Σ⁺ states of CO molecule in aug-cc-pVDZ basis using a bond distance of 1.129 Å. Corrections at a given order, ω^(k), and the sum of excitation energies for both states at a given order, ω_k(F¹Σ⁺) + ω_k(4¹Σ⁺), are also reported. Note that ω⁽⁰⁾ = ω^CCS, ω⁽¹⁾ = 0, and ω₂ = ω^CIS(D). All results in eV.

	F¹Σ⁺		4¹Σ⁺
k	ω^(k)	ω_k	ω^(k)	ω_k	ω_k(F¹Σ⁺) + ω_k(4¹Σ⁺)
0	15.967	15.967	16.524	16.524	32.491
1	0.000	15.967	0.000	16.524	32.491
2	−1.416	14.551	−1.745	14.779	29.330
3	0.207	14.758	0.309	15.089	29.846
4	0.054	14.812	−0.083	15.006	29.817
5	0.130	14.942	0.069	15.075	30.016
6	−0.258	14.683	0.178	15.253	29.936
7	0.120	14.803	−0.158	15.095	29.898
8	−0.215	14.588	0.251	15.346	29.933
9	0.201	14.788	−0.239	15.106	29.895
10	−0.195	14.593	0.231	15.338	29.931
11	0.136	14.729	−0.151	15.187	29.916
12	−0.051	14.678	0.052	15.239	29.917
13	−0.039	14.640	0.051	15.290	29.930
14	0.199	14.839	−0.228	15.062	29.901
15	−0.337	14.502	0.377	15.439	29.941
16	0.527	15.028	−0.575	14.864	29.892
17	−0.656	14.372	0.710	15.574	29.946
18	0.777	15.149	−0.832	14.742	29.891
19	−0.792	14.357	0.843	15.585	29.942
20	0.681	15.038	−0.724	14.861	29.898
21	−0.401	14.637	0.434	15.295	29.932
22	−0.115	14.522	0.094	15.389	29.911
23	0.837	15.359	−0.830	14.559	29.918
24	−1.818	13.541	1.823	16.383	29.923
25	2.925	16.466	−2.941	13.442	29.908
26	−4.091	12.375	4.114	17.556	29.930
27	5.057	17.432	−5.083	12.473	29.905
28	−5.558	11.874	5.583	18.056	29.930
29	5.218	17.091	−5.238	12.818	29.910
30	−3.603	13.488	3.615	16.434	29.922
31	0.356	13.844	−0.358	16.076	29.920
32	4.875	18.719	−4.884	11.192	29.911
33	−12.089	6.630	12.109	23.300	29.930
34	21.025	27.655	−21.054	2.246	29.901
35	−30.757	−3.102	30.793	33.039	29.937
36	39.740	36.638	−39.780	−6.740	29.898
37	−45.579	−8.940	45.618	38.877	29.937
38	45.087	36.147	−45.122	−6.245	29.902
39	−34.534	1.613	34.561	28.316	29.929
40	9.951	11.563	−9.967	18.349	29.912
41	31.900	43.463	−31.895	−13.546	29.917

	F¹Σ⁺		4¹Σ⁺
k	ω^(k)	ω_k	ω^(k)	ω_k	ω_k(F¹Σ⁺) + ω_k(4¹Σ⁺)
0	15.967	15.967	16.524	16.524	32.491
1	0.000	15.967	0.000	16.524	32.491
2	−1.416	14.551	−1.745	14.779	29.330
3	0.207	14.758	0.309	15.089	29.846
4	0.054	14.812	−0.083	15.006	29.817
5	0.130	14.942	0.069	15.075	30.016
6	−0.258	14.683	0.178	15.253	29.936
7	0.120	14.803	−0.158	15.095	29.898
8	−0.215	14.588	0.251	15.346	29.933
9	0.201	14.788	−0.239	15.106	29.895
10	−0.195	14.593	0.231	15.338	29.931
11	0.136	14.729	−0.151	15.187	29.916
12	−0.051	14.678	0.052	15.239	29.917
13	−0.039	14.640	0.051	15.290	29.930
14	0.199	14.839	−0.228	15.062	29.901
15	−0.337	14.502	0.377	15.439	29.941
16	0.527	15.028	−0.575	14.864	29.892
17	−0.656	14.372	0.710	15.574	29.946
18	0.777	15.149	−0.832	14.742	29.891
19	−0.792	14.357	0.843	15.585	29.942
20	0.681	15.038	−0.724	14.861	29.898
21	−0.401	14.637	0.434	15.295	29.932
22	−0.115	14.522	0.094	15.389	29.911
23	0.837	15.359	−0.830	14.559	29.918
24	−1.818	13.541	1.823	16.383	29.923
25	2.925	16.466	−2.941	13.442	29.908
26	−4.091	12.375	4.114	17.556	29.930
27	5.057	17.432	−5.083	12.473	29.905
28	−5.558	11.874	5.583	18.056	29.930
29	5.218	17.091	−5.238	12.818	29.910
30	−3.603	13.488	3.615	16.434	29.922
31	0.356	13.844	−0.358	16.076	29.920
32	4.875	18.719	−4.884	11.192	29.911
33	−12.089	6.630	12.109	23.300	29.930
34	21.025	27.655	−21.054	2.246	29.901
35	−30.757	−3.102	30.793	33.039	29.937
36	39.740	36.638	−39.780	−6.740	29.898
37	−45.579	−8.940	45.618	38.877	29.937
38	45.087	36.147	−45.122	−6.245	29.902
39	−34.534	1.613	34.561	28.316	29.929
40	9.951	11.563	−9.967	18.349	29.912
41	31.900	43.463	−31.895	−13.546	29.917

FIG. 12.

View large Download slide

Excitation energy corrections, ω^(k) (eV), for the F¹Σ⁺ and 4¹Σ⁺ states of the CO molecule calculated using CPS(D-n) series through order 40 in aug-cc-pVDZ basis using a bond distance of 1.129 Å. The horizontal line, where the excitation energy is zero, is a mirror plane for the higher-order excitation energy corrections.

The sum of the CCSD excitation energy for the F¹Σ⁺ and the 4¹Σ⁺ state is

ω_{F}^{C C S D} + ω_{4}^{C C S D} = 29.92 eV .

This sum is clearly obtained at about order 10 in the last column of Table III and maintained at the consecutive orders, despite the excitation energy corrections becoming larger than the total excitation energies and the calculated excitation energies at some orders becoming negative.

For the (B¹Σ⁺, C¹Σ⁺) and (F¹Σ⁺, 4¹Σ⁺) excitation pairs, a front- and a back-door intruder state, respectively, gave rise to the divergent pairs of ripple pattern. For the perturbation series reported in Paper II,¹⁸ we have also seen examples of excitation pairs where the intruder state for the excitation pair is replaced with a crosser state and where the overlaying divergent ripple patterns for the excitation pair therefore become replaced with overlaying convergent ripple patterns.

The above development shows that the asymptotic convergence of the excitation energy series can be reproduced by the two-state model. In fact, if the convergence/divergence behavior of the excitation energy series is examined in detail and through high enough order, we expect that the asymptotic convergence patterns and the convergence rate of the series will be dictated by a two-state problem, except for pathological cases.

VI. CONVERGENCE OF CP SERIES FOR MOLECULAR PROPERTIES

We now consider the convergence of CP series for molecular properties. In Subsection VI A, we consider the convergence of the CP series for first order molecular properties, and in Subsection VI B, the convergence is considered for molecular properties that can be described by the linear response function. Subsection VI C gives examples of the convergence of CP series for the various molecular properties.

A. First-order molecular properties

The first-order molecular property for an operator X and a CC target state can be written as¹⁷

{⟨ ⟨ X ⟩ ⟩}_{0} = ⟨) H F | X^{T} | H F (⟩ + \sum_{i = 1}^{t} \sum_{μ_{i}} {\bar{t}}_{μ_{i}} ⟨) μ_{i} | X^{T} | H F (⟩ .

(89)

The cluster amplitudes satisfy the cluster amplitude equations,

⟨) μ_{i} | H_{0}^{T} | H F (⟩ = 0, i = 1,2, \dots, t,

(90)

and the multipliers ${\bar{t}}_{μ_{i}}$ satisfy the multiplier equation,

\bar{t} J = - η,

(91)

where the right-hand side is given by

η_{μ_{i}} = ⟨) H F | H_{0}^{T} | μ_{i} (⟩,

(92)

and the elements of the Jacobian matrix for the CC target state, J, read

J_{μ_{i} ν_{j}} = ⟨) μ_{i} | [) H_{0}^{T}, θ_{ν_{j}}]) | H F (⟩, i, j = 1,2, \dots, t .

(93)

In CP theory, the CC target state is parameterized with the CC parent state as an expansion point. This implies that the similarity-transformed operators O^T in Eqs. (89)–(93) are replaced by $e^{- δ T} O^{^{*} T} e^{δ T}$ ⁠. For example, the Jacobian in Eq. (93) becomes the Jacobian in Eq. (45), and the amplitude equations in Eq. (90) become the amplitude equations in Eqs. (41) and (42). The cluster amplitudes δt and the multipliers $\bar{t}$ can be expanded in orders of $Φ^{^{*} T}$ ⁠, and CP series can be determined for these parameters. The CP series for δt is given in Eq. (55), and the CP series for the multipliers becomes

{\bar{t}}_{μ_{i}} = {\bar{t}}_{μ_{i}}^{(1)} + {\bar{t}}_{μ_{i}}^{(2)} + \dots,

(94)

with a vanishing zeroth-order term. Substituting these expansions into the first-order molecular property expression in Eq. (89), we obtain a CP series for the first-order molecular property,

{⟨ ⟨ X ⟩ ⟩}_{0} = {⟨ ⟨ X ⟩ ⟩}_{0}^{P} + {⟨ ⟨ X ⟩ ⟩}_{0}^{(1)} + {⟨ ⟨ X ⟩ ⟩}_{0}^{(2)} + \dots,

(95)

where ${⟨ ⟨ X ⟩ ⟩}_{0}^{P}$ is the first-order molecular property for the CC parent state and ${⟨ ⟨ X ⟩ ⟩}_{0}^{(k)}$ is the kth-order term.

We now consider the fictitious system where the complex strength parameter z is multiplied on the perturbation $Φ^{^{*} T}$ ⁠. The first-order molecular property ⟨⟨X⟩⟩₀(z) becomes an analytic function of z, δt(z), and $\bar{t} (z)$ in the complex plane. For z = 0, ⟨⟨X⟩⟩₀(z) becomes the first-order molecular property for the CC parent state ${⟨ ⟨ X ⟩ ⟩}_{0}^{P}$ ⁠, and for z = 1, ⟨⟨X⟩⟩₀(z) becomes the first-order molecular property for the physical system given in Eq. (89). ⟨⟨X⟩⟩₀(z) can be expanded with z = 0 as expansion point, giving

{⟨ ⟨ X ⟩ ⟩}_{0} (z) = {⟨ ⟨ X ⟩ ⟩}_{0}^{P} + \sum_{n = 1}^{\infty} z^{n} () () X)) {))}_{0}^{(n)} .

(96)

For z = 1, Eq. (96) gives the CP series in Eq. (95). For the CP series to be convergent, we must require that the CP series for the cluster amplitudes converge and therefore that the Jacobian J(z) has no singularities in the complex plane inside the unit circle |z| ≤ 1. Furthermore, the first-order molecular property ⟨⟨X⟩⟩₀(z) for the fictitious system must not have any singularities inside the unit circle. The singularities of ⟨⟨X⟩⟩₀(z) enter through the singularities in the z dependence of the multipliers $\bar{t} (z)$ which in turn require that J(z) has no singularities inside the unit circle |z| ≤ 1. This condition is equivalent to the condition for having a convergent CP series for the cluster amplitudes.

When ${⟨ ⟨ X ⟩ ⟩}_{0}^{(1)}$ in Eq. (95) is evaluated, it contains terms that are linear in δt⁽¹⁾ and ${\bar{t}}^{(1)}$ ⁠, whereas the terms that are bi-linear in δt⁽¹⁾ and ${\bar{t}}^{(1)}$ enter in second order. When evaluating ${⟨ ⟨ X ⟩ ⟩}_{0}^{(1)}$ ⁠, the leading-order computational scaling arises from the determination of δt⁽¹⁾ and ${\bar{t}}^{(1)}$ ⁠. The terms that are bi-linear in δt⁽¹⁾ and ${\bar{t}}^{(1)}$ can be calculated without any extra cost when the ${⟨ ⟨ X ⟩ ⟩}_{0}^{(1)}$ is evaluated, and it therefore becomes attractive to restructure the series in Eq. (95) by introducing a generalized order concept where the term bi-linear in δt⁽¹⁾ and ${\bar{t}}^{(1)}$ is accounted for as a first-order contribution. To implement this generalized order concept for an arbitrary order k, we introduce the total kth-order amplitudes and multipliers,

^{k} δ t = δ t^{(1)} + δ t^{(2)} + \dots + δ t^{(k)},

(97)

^{k} \bar{t} = {\bar{t}}^{(1)} + {\bar{t}}^{(2)} + \dots + {\bar{t}}^{(k)},

(98)

and evaluate ⟨⟨X⟩⟩₀through kth order using these total amplitudes and multipliers,

^{k} {⟨ ⟨ X ⟩ ⟩}_{0} = ⟨) H F |^{k} X^{^{*} T} | H F (⟩ + \sum_{i = 1}^{t} \sum_{μ_{i}}^{k} {\bar{t}}_{μ_{i}} ⟨) μ_{i} |^{k} X^{^{*} T} | H F (⟩,

(99)

where

^{k} X^{^{*} T} = X^{^{*} T} + [) X^{^{*} T},^{k} δ T]) + \frac{1}{2} [) [) X^{^{*} T},^{k} δ T]),^{k} δ T])

(100)

and

^{k} δ T = \sum_{μ_{i}}^{k} {δ t}_{μ_{i}} θ_{μ_{i}} .

(101)

Note that, formally,

lim_{k \to \infty}^{k} {⟨ ⟨ X ⟩ ⟩}_{0} = {⟨ ⟨ X ⟩ ⟩}_{0} .

(102)

The leading-order computational scaling for evaluating ^k⟨⟨X⟩⟩₀ arises from the determination of δt^(k) and ${\bar{t}}^{(k)}$ ⁠, as for ${⟨ ⟨ X ⟩ ⟩}_{0}^{(k)}$ ⁠. The major difference between evaluating ⟨⟨X⟩⟩₀ using the series in Eq. (95) and the one defined by Eq. (99) is that using Eq. (95), contributions are evaluated strictly of order k. By contrast, when Eq. (99) is used, the total contribution through order k is evaluated, but Eq. (99) in addition also contains all higher-order contributions that can be generated using the amplitude and multiplier corrections through order k. When the series in Eq. (95) is convergent, the series defined by Eq. (99) will also be convergent as the series in Eq. (99) results from a simple restructuring of the terms in the series in Eq. (95).

B. The linear response function

For the operator-frequency pairs (X, ω_X) and (Y, ω_Y), the linear response function for a CC target state can be written as¹⁷

\begin{matrix} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} = \frac{1}{2} C^{\pm ω_{Y}} P_{ω_{X} ω_{Y}}^{X Y} () η^{X} t^{Y} (ω_{Y}) + \frac{1}{2} {() t^{X} (ω_{X})))}^{T} F t^{Y} (ω_{Y}))), \\ ω_{X} + ω_{Y} = 0, \end{matrix}

(103)

where

η_{μ_{i}}^{X} = ⟨) H F | X^{T} | μ_{i} (⟩ + \sum_{j} \sum_{ν_{j}} {\bar{t}}_{ν_{j}} ⟨) ν_{j} | [) X^{T}, θ_{μ_{i}}]) | H F (⟩,

(104)

\begin{matrix} F_{μ_{i} ν_{j}} & = & ⟨) H F | [) [) Φ, θ_{μ_{i}}]), θ_{ν_{j}}]) | H F (⟩ δ_{i 1} δ_{j 1} \\ + \sum_{m} \sum_{λ_{m}} {\bar{t}}_{λ_{m}} ⟨) λ_{m} | [) [) Φ^{T}, θ_{μ_{i}}]), θ_{ν_{j}}]) | H F (⟩ . \end{matrix}

(105)

$P_{ω_{X} ω_{Y}}^{X Y}$ is the (operator, frequency) pair permutation operator, and C^±ω is the frequency permutation operator. The cluster amplitudes satisfy the cluster amplitude equation in Eq. (90), and the multipliers satisfy Eq. (91). The first-order frequency dependent amplitude responses t^Y(ω_Y) are determined from the response amplitude equations,

() J - ω_{X} I)) t^{X} (ω_{X}) = - ξ^{X},

(106)

where the right-hand side is given by

ξ_{μ_{i}}^{X} = ⟨) μ_{i} | X^{T} | H F (⟩ .

(107)

In CP theory, the CC target state is parameterized with the CC parent state as an expansion point and the similarity-transformed operators O^T in Eqs. (103)–(107) therefore have to be replaced by $e^{- δ T} O^{^{*} T} e^{δ T}$ as described in Subsection VI A. The cluster amplitudes δt, the multipliers $\bar{t}$ ⁠, and the response amplitudes t^X(ω_X) can be expanded in orders of $Φ^{^{*} T}$ ⁠, and CP series can be determined for these parameters. The CP series for δt is given in Eq. (55) and for $\bar{t}$ in Eq. (94). The CP perturbation series for t^X(ω_X) is derived in Ref. 28. Substituting these expansions into the linear response function, we can determine a CP series for the linear response function,

{⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} = {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{P} + \sum_{n = 1}^{\infty} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{(n)},

(108)

where ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{P}$ is the linear response function for the CC parent state.²⁸

We now consider the fictitious system where the complex strength parameter z is multiplied on the perturbation $Φ^{^{*} T}$ ⁠. The linear response function for the fictitious system ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} (z)$ becomes an analytic function of z in the complex plane that depends on δt, $\bar{t}$ and t^X(ω_X) and can be expanded with z = 0 as an expansion point, giving

{⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} (z) = {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{P} + \sum_{n = 1}^{\infty} z^{n} () () X; Y)) {))}_{ω_{Y}}^{(n)} (z) .

(109)

For z = 1, Eq. (109) gives the CP series for the linear response function in Eq. (108). For the CP series in Eq. (108) to converge, the CP series for the cluster amplitudes must converge implying that the Jacobian J(z) has no singularity in the complex plane inside the unit circle |z| ≤ 1. Furthermore, ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} (z)$ must not have singularities in the complex plane inside the unit circle |z| ≤ 1. The singularities of ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}} (z)$ enter through singularities in z dependence of the multipliers $\bar{t}$ and the response amplitudes t^Y(ω_Y). For the multipliers, we must require that J(z) has no singularities inside the unit circle |z| ≤ 1. This condition is equivalent to the condition for having a convergent series for the cluster amplitudes. For the response amplitudes, we must also require that J(z) − ω_YI has no singularities inside the unit circle |z| ≤ 1. This condition restricts the number of convergent expansions of the linear response function, in particular when ω_Y is close to an eigenvalue ω_x of the Jacobian J.

When ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{(k)}$ is evaluated, it requires the amplitudes and multipliers, δt^(k), ${\bar{t}}^{(k)}$ ⁠, and ${t^{X}}^{(k)} (ω_{X})$ ⁠. The leading-order computational scaling for evaluating ${⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}^{(k)}$ arises from the evaluation of δt^(k), ${\bar{t}}^{(k)}$ and ${t^{X}}^{(k)} (ω_{X})$ ⁠. As for the evaluation of ${⟨ ⟨ X ⟩ ⟩}_{0}^{(k)}$ ⁠, it becomes attractive to introduce a generalized order concept, where the order k embraces all the terms that can be determined from the amplitude, multiplier, and response vector corrections, δt^(p), ${\bar{t}}^{(p)}$ ⁠, and ${t^{X}}^{(p)} (ω_{X})$ ⁠, where p = 0, 1, …, k. We denote the kth-order term in this series as $^{k} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}$ ⁠. In practice, the evaluation of this kth-order term requires that we introduce

^{k} t^{X} (ω_{X}) = {t^{X}}^{(1)} (ω_{X}) + {t^{X}}^{(2)} (ω_{X}) + \dots + {t^{X}}^{(k)} (ω_{X})

(110)

in addition to ^kδt and $^{k} \bar{t}$ of Eqs. (97) and (98), respectively. Using these total kth-order amplitudes and multipliers, $^{k} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}$ can be evaluated from Eq. (103) in the same way ^k⟨⟨X⟩⟩₀ was evaluated from Eq. (89) in Eq. (99). In Ref. 28, details are given on how $^{k} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}$ is determined. When the series in Eq. (108) is convergent, the $^{k} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}$ series will also be convergent as the $^{k} {⟨ ⟨ X; Y ⟩ ⟩}_{ω_{Y}}$ series results from a simple restructuring of the terms in the series in Eq. (108).

C. Numerical illustrations of convergence of CP series for molecular properties

As an illustration of the convergence requirements for CP series for molecular properties, we consider calculations on CO for R_CO = 1.129 Å, using the aug-cc-pVDZ basis. The CPS(D-n) polarizability and dipole moment corrections through an arbitrary order have been implemented using Psithon—a Python interface to the PSI4 program²⁴—and employing the Numerical Python Library (Numpy).²⁵ The convergence of the CPS(D-n) series for the Σ component of the static polarizability calculated using Eq. (108), the dipole moment calculated using Eq. (95), and the energy calculated as described in Paper I¹ are displayed in Fig. 13. The asymptotic convergence of the static polarizability, the dipole moment, and the energy show a geometric progression. The convergence rate of the CP series for the static polarizability, the dipole moment, and the energy are very similar. For the energy, the convergence rate is determined by the inverse distance from the unit circle for the crosser state associated with the primary critical point and this appears also to be the case for the static polarizability and the dipole moment.

FIG. 13.

View large Download slide

Absolute CPS(D) errors (in atomic units), as compared to target CCSD values, for the ground-state energy, E₀, Σ component of the permanent dipole moment, ⟨Σ⟩₀, and Σ component of the static dipole polarizability, ⟨⟨Σ; Σ⟩⟩₀, calculated through a given order for the CO molecule in aug-cc-pVDZ basis using a bond distance of 1.129 Å.

VII. SUMMARY AND CONCLUSION

In CP theory, perturbation series are determined in orders of the CC parent state similarity transformed fluctuation potential where the zeroth order term in the series is the energy or molecular property for the CC parent state and where the perturbation series formally converge to the energy or molecular property for the CC target state. In this paper, we have investigated the theoretical foundation for having convergent CP series for the ground-state energy and for molecular properties including excitation energies.

The theoretical foundation for convergent CP series is established by considering a fictitious system where a complex scaling parameter z is multiplied on the perturbation. The condition for having a convergent CP ground-state energy series is that the derivative of the amplitude equation for the fictitious system is non-singular inside the unit circle |z| ≤ 1. This condition is equivalent to requiring that the Jacobian J(z) for the fictitious system is non-singular within the unit circle |z| ≤ 1.

To have a convergent CP series for an excitation energy ω_x, we must in addition require that the excitation energy for the fictitious system, ω_x(z), does not become degenerate in the complex plane inside the unit circle |z| ≤ 1. Since ω_x(z) is determined by diagonalizing the Jacobian J(z) for the fictitious system, the condition for convergence becomes that there must not be another excitation energy ω_y(z) that is equal to ω_x(z) inside the unit circle |z| ≤ 1.

To have a convergent CP molecular property series, the molecular property expression for the fictitious system must not have any singularities inside the unit circle |z| ≤ 1. Singularities in molecular property expressions are caused by singularities in the cluster amplitudes, multipliers, or response vectors. For the cluster amplitudes and multipliers, singularities can arise in the Jacobian J(z). For frequency dependent molecular properties, additional singularities can arise through singularities in the frequency-shifted Jacobian.

The identification of singularities and degeneracies in the Jacobian for the fictitious system can be performed by searching for avoided crossings for the Jacobian J(z) on the real axis $Re (z)$ ⁠. For the ground-state energy, the search for singularities identifies mechanisms giving front-door and back-door intruder states. The two states composing the avoided crossing for the singularity with the smallest distance to the origin can be used to set up a two-state model for the fictitious system, from which we can obtain an accurate description of the convergence rate and the convergence patterns of the higher-order terms in the CP ground-state energy series. Numerical examples are given, where the convergence rate and the convergence patterns of the higher-order terms in the CP ground-state energy series are modeled using the two-state model.

When examining the convergence for an excitation energy ω_x, we look for degeneracies and thus avoided crossings between ω_x(z) and another excitation energy ω_y(z), when the Jacobian response eigenvalue equation is solved on the real axis $Re (z)$ ⁠. For the two excitations that compose the avoided crossing with the smallest distance to the origin, we can set up a two-state excitation model for the Jacobian eigenvalue equation for the fictitious system, and from this two-state excitation model, we can obtain an accurate description of the convergence rate and the convergence patterns of the higher-order terms in the CP excitation energy series. We give numerical examples where the convergence patterns of CP excitation energy series is modeled using the two-state excitation model. We also consider a CP excitation energy pair, ω_x(z) and ω_y(z), where ω_y(z) is an intruder state for ω_x(z) and vice versa, and where the CP series for these two excitation energies therefore diverge. For such an excitation energy pair, we have shown that the excitation energy corrections for higher orders have opposite sign, and further that the sum of the excitation energies for higher orders becomes equal to the sum of the target excitation energies, provided that the two-state model accurately describes the higher-order terms in the CP excitation energy series. We have reported numerical calculations, which show that the above relations are fulfilled and that the two-state model therefore gives an accurate description of the higher-order terms in the CP excitation energy series. We also note that the point of degeneracies ω_x(z^c) = ω_y(z^c) inside the unit circle for ω_x(z) and ω_y(z) does not directly affect the convergence of the other excitation energies. Summarizing, the calculations of CP series for the ground-state energy and for excitation energies show that the convergence of CP series for the ground-state energy and for excitation energies effectively becomes the convergence of a two-state problem at higher orders. We expect this to be a general behavior, except for pathological cases.

ACKNOWLEDGMENTS

J.O. acknowledges support from the Danish Council for Independent Research, Grant No. DFF-4181-00537. This research used resources of the Oak Ridge Leadership Computing Facility at the Oak Ridge National Laboratory, which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC05-00OR22725 and resources provided by the Wroclaw Centre for Networking and Supercomputing (http://www.wcss.pl/en/), Grant No. 277.

REFERENCES

1.

F.

Pawłowski

,

J.

Olsen

, and

P.

Jørgensen

, “

Cluster perturbation theory. I. Theoretical foundation for a coupled cluster target state and ground-state energies

,”

J. Chem. Phys.

150

,

134108

(

2019

).

https://doi.org/10.1063/1.5004037

Google Scholar

Crossref

PubMed

2.

C.

Møller

and

M. S.

Plesset

,

Phys. Rev.

46

,

618

(

1934

).

https://doi.org/10.1103/physrev.46.618

Google Scholar

Crossref

3.

D.

Cremer

,

Wiley Interdiscip. Rev.: Comput. Mol. Sci.

1

,

509

(

2011

).

https://doi.org/10.1002/wcms.58

Google Scholar

Crossref

4.

J. J.

Eriksen

,

K.

Kristensen

,

T.

Kjærgaard

,

P.

Jørgensen

, and

J.

Gauss

,

J. Chem. Phys.

140

,

064108

(

2014

).

https://doi.org/10.1063/1.4862501

Google Scholar

Crossref

PubMed

5.

J. J.

Eriksen

,

P.

Jørgensen

, and

J.

Gauss

,

J. Chem. Phys.

142

,

014102

(

2015

).

https://doi.org/10.1063/1.4904754

Google Scholar

Crossref

PubMed

6.

K.

Kristensen

,

J. J.

Eriksen

,

D. A.

Matthews

,

J.

Olsen

, and

P.

Jørgensen

,

J. Chem. Phys.

144

,

064103

(

2016

).

https://doi.org/10.1063/1.4941605

Google Scholar

Crossref

PubMed

7.

J.

Olsen

,

P.

Jørgensen

,

T.

Helgaker

, and

O.

Christiansen

,

J. Chem. Phys.

112

,

9736

(

2000

).

https://doi.org/10.1063/1.481611

Google Scholar

Crossref

8.

T.

Helgaker

,

P.

Jørgensen

, and

J.

Olsen

,

Molecular Electronic-Structure Theory

(

Wiley

,

Chichester

,

2000

).

Google Scholar

Crossref

9.

J.

Olsen

and

P.

Jørgensen

, “

Convergence patterns and rates in two-state perturbation expansions

,”

J. Chem. Phys.

(submitted).

10.

J.

Olsen

,

O.

Christiansen

,

H.

Koch

, and

P.

Jørgensen

,

J. Chem. Phys.

105

,

5082

(

1996

).

https://doi.org/10.1063/1.472352

Google Scholar

Crossref

11.

O.

Christiansen

,

J.

Olsen

,

P.

Jørgensen

,

H.

Koch

, and

P.-Å.

Malmqvist

,

Chem. Phys. Lett.

261

,

369

(

1996

).

https://doi.org/10.1016/0009-2614(96)00974-8

Google Scholar

Crossref

12.

M. L.

Leininger

,

W. D.

Allen

,

H. F.

Schaefer

, and

C. D.

Sherrill

,

J. Chem. Phys.

112

,

9213

(

2000

).

https://doi.org/10.1063/1.481764

Google Scholar

Crossref

13.

J. J.

Eriksen

,

K.

Kristensen

,

D. A.

Matthews

,

P.

Jørgensen

, and

J.

Olsen

,

J. Chem. Phys.

145

,

224104

(

2016

).

https://doi.org/10.1063/1.4971294

Google Scholar

Crossref

PubMed

14.

N. C.

Handy

,

P. J.

Knowles

, and

K.

Somasundram

,

Theor. Chem. Acc.

68

,

87

(

1985

).

https://doi.org/10.1007/bf00698753

Google Scholar

Crossref

15.

J. F.

Stanton

and

R. J.

Bartlett

,

J. Chem. Phys.

98

,

7029

(

1993

).

https://doi.org/10.1063/1.464746

Google Scholar

Crossref

16.

I.

Shavitt

and

R. J.

Bartlett

,

Many-Body Methods in Chemistry and Physics: MBPT and Coupled-Cluster Theory

(

Cambridge University Press

,

2009

).

Google Scholar

Crossref

17.

H.

Koch

and

P.

Jørgensen

,

J. Chem. Phys.

93

,

3333

(

1990

).

https://doi.org/10.1063/1.458814

Google Scholar

Crossref

18.

F.

Pawłowski

,

J.

Olsen

, and

P.

Jørgensen

, “

Cluster perturbation theory. II. Excitation energies for a coupled cluster target state

,”

J. Chem. Phys.

150

,

134109

(

2019

).

https://doi.org/10.1063/1.5053167

Google Scholar

Crossref

PubMed

19.

K.

Fritzsche

and

H.

Grauert

,

From Holomorphic Functions to Complex Manifolds, Graduate Texts in Mathematics

(

Springer-Verlag

,

New York

,

2002

), p.

34

.

Google Scholar

Crossref

20.

T. H.

Dunning

, Jr.,

J. Chem. Phys.

90

,

1007

(

1989

).

https://doi.org/10.1063/1.456153

Google Scholar

Crossref

21.

R. A.

Kendall

,

T. H.

Dunning

, Jr., and

R. J.

Harrison

,

J. Chem. Phys.

96

,

6796

(

1992

).

https://doi.org/10.1063/1.462569

Google Scholar

Crossref

22.

J.

Olsen

,

Lucia—A Configuration Interaction and Coupled Cluster Program

(

University of Aarhus

,

2016

).

Google Scholar

23.

J.

Olsen

,

J. Chem. Phys.

113

,

7140

(

2000

).

https://doi.org/10.1063/1.1290005

Google Scholar

Crossref

24.

R. M.

Parrish

,

L. A.

Burns

,

D. G.

Smith

,

A. C.

Simmonett

,

A. E.

DePrince

,

E. G.

Hohenstein

,

U.

Bozkaya

,

A. Y.

Sokolov

,

R.

Di Remigio

,

R. M.

Richard

,

J. F.

Gonthier

,

A. M.

James

,

H. R.

McAlexander

,

A.

Kumar

,

M.

Saitow

,

X.

Wang

,

B. P.

Pritchard

,

P.

Verma

,

H. F.

Schaefer

,

K.

Patkowski

,

R. A.

King

,

E. F.

Valeev

,

F. A.

Evangelista

,

J. M.

Turney

,

T. D.

Crawford

, and

C. D.

Sherrill

,

J. Chem. Theory Comput.

13

,

3185

(

2017

).

https://doi.org/10.1021/acs.jctc.7b00174

Google Scholar

Crossref

PubMed

25.

S.

Van Der Walt

,

S. C.

Colbert

, and

G.

Varoquaux

,

Comput. Sci. Eng.

13

,

22

(

2011

); e-print arXiv:1102.1523.

https://doi.org/10.1109/mcse.2011.37

26.

P.

Baudin

,

F.

Pawłowski

,

D.

Bykov

,

D.

Liakh

,

K.

Kristensen

,

J.

Olsen

, and

P.

Jørgensen

, “

Cluster perturbation theory. III. Perturbation series for coupled cluster singles and doubles excitation energies

,”

J. Chem. Phys.

150

,

134110

(

2019

).

https://doi.org/10.1063/1.5046935

Google Scholar

Crossref

PubMed

27.

M.

Head-Gordon

,

R. J.

Rico

,

M.

Oumi

, and

T. J.

Lee

,

Chem. Phys. Lett.

219

,

21

(

1994

).

https://doi.org/10.1016/0009-2614(94)00070-0

Google Scholar

Crossref

28.

F.

Pawłowski

,

D.

Bykov

,

D.

Liakh

,

J.

Olsen

, and

P.

Jørgensen

, “

Cluster perturbation theory. VI. Linear response function for a coupled cluster target state

” (unpublished).

2019

Author(s)

Cluster perturbation theory. IV. Convergence of cluster perturbation series for energies and molecular properties

I. INTRODUCTION

II. PERTURBATION THEORY WITH A MØLLER-PLESSET PARTITIONING OF THE HAMILTONIAN

A. Requirements for convergent energy perturbation series

B. The two-state model

1. The two-state model and its critical points

2. Perturbation expansions of the energies for symmetric perturbations

3. Archetypes of convergence

C. Outline

III. CLUSTER PERTURBATION THEORY FOR THE ENERGY AND AMPLITUDE EQUATIONS AND FOR THE JACOBIAN

IV. CONVERGENCE OF THE CP ENERGY SERIES

A. Theory

B. Numerical illustrations

V. CONVERGENCE OF CP EXCITATION ENERGY SERIES

A. Theory

B. Numerical illustrations

VI. CONVERGENCE OF CP SERIES FOR MOLECULAR PROPERTIES

A. First-order molecular properties

B. The linear response function

C. Numerical illustrations of convergence of CP series for molecular properties

VII. SUMMARY AND CONCLUSION

ACKNOWLEDGMENTS

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

Cluster perturbation theory. IV. Convergence of cluster perturbation series for energies and molecular properties

I. INTRODUCTION

II. PERTURBATION THEORY WITH A MØLLER-PLESSET PARTITIONING OF THE HAMILTONIAN

A. Requirements for convergent energy perturbation series

B. The two-state model

1. The two-state model and its critical points

2. Perturbation expansions of the energies for symmetric perturbations

3. Archetypes of convergence

C. Outline

III. CLUSTER PERTURBATION THEORY FOR THE ENERGY AND AMPLITUDE EQUATIONS AND FOR THE JACOBIAN

IV. CONVERGENCE OF THE CP ENERGY SERIES

A. Theory

B. Numerical illustrations

V. CONVERGENCE OF CP EXCITATION ENERGY SERIES

A. Theory

B. Numerical illustrations

VI. CONVERGENCE OF CP SERIES FOR MOLECULAR PROPERTIES

A. First-order molecular properties

B. The linear response function

C. Numerical illustrations of convergence of CP series for molecular properties

VII. SUMMARY AND CONCLUSION

ACKNOWLEDGMENTS

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only