Guaranteed convergence for a class of coupled-cluster methods based on Arponen's extended theory

Abstract

A wide class of coupled-cluster methods is introduced, based on Arponen's extended coupled-cluster theory. This class of methods is formulated in terms of a coordinate transformation of the cluster operators. The mathematical framework for the error analysis of coupled-cluster methods based on Arponen's bivariational principle is presented, in which the concept of local strong monotonicity of the flipped gradient of the energy is central. A general mathematical result is presented, describing sufficient conditions for coordinate transformations to preserve the local strong monotonicity. The result is applied to the presented class of methods, which include the standard and quadratic coupled-cluster methods, and also Arponen's canonical version of extended coupled-cluster theory. Some numerical experiments are presented, and the use of canonical coordinates for diagnostics is discussed.

GRAPHICAL ABSTRACT

Keywords:

• Coupled-cluster method
• extended coupled-cluster method
• error analysis
• electronic-structure theory

1. Introduction

It is with delight that the authors dedicate this work to Professor Jürgen Gauß on the occasion of his sixtieth birthday. In the spirit of his pursuit of scientific rigour, especially the attention to detail in coupled-cluster (CC) theory, we here present a mathematical study of some alternative formulations based on Arponen's extended CC (ECC) method [1,2]. The ECC method is defined in terms of the critical points of an energy functional, (1) ${\mathcal{E}}_{\mathrm{ECC}}(T,\mathrm{\Lambda })=⟨{\varphi }_0,e^{{\mathrm{\Lambda }}^{†}}{e}^{-T}H{e}^T{\varphi }_0⟩,$(1) where T and Λ are cluster operators in the usual sense of CC theory. The well-known CC Lagrangian introduced by Helgaker and Jørgensen [3] is obtained by expanding $e^{{\mathrm{\Lambda }}^{†}}$ to first order. In this sense, the standard CC approach is an approximation to the ECC method.

We will study a collection of methods that generalises this idea, defined by substitution of $e^{{\mathrm{\Lambda }}^{†}}$ by a Taylor polynomial of fixed degree n. This can be formulated as a coordinate transformation of the cluster amplitudes. A second class of models is obtained by a further coordinate transformation introduced by Arponen to ensure a certain linkedness structure of the energy functional. We refer to these coordinates as canonical, as the time-dependent Schrödinger equation takes the form of Hamilton's equations of motion in this case. Thus, we obtain two hierarchies NC-ECC$(n)$, using non-canonical coordinates, and C-ECC$(n)$ using canonical coordinates. Our mathematical results imply that when the cluster operators are not truncated, all these models are exact and equivalent to the time-independent Schrödinger equation. Moreover, Galerkin approximations (i.e. generic truncation schemes that can approach the untruncated limit) will converge under certain relatively mild single-reference-type conditions. While the methods discussed here are all expensive (for n>1), we consider them as stepping stones towards producing competitive alternatives to standard CC theory that alleviate deficiencies of the latter, such as the inability to correctly break chemical bonds.

The various forms of CC methods are today among the most widely used for wavefunction-based calculations on manybody systems. The main idea stems from Hubbard's exponential parameterisation of the wavefunction based on cluster operators in manybody perturbation theory [4], which was taken as starting point for ab initio treatments by Coester and Kümmel for nuclear structure calculations in the 1950s [5,6]. The modern form of standard CC theory was developed by, among others, Sinanoğlu, Paldus and Cizek in the 1960s [7] and the CC method with singles, doubles and perturbative triples [CCSD(T)] today constitutes ‘the gold standard of quantum chemistry’ due to its excellent balance between computational cost and accuracy [8]. In nuclear structure calculations the same method has gained traction in the last decade, providing excellent predictive power for light to medium nuclei [9]. Coupled-cluster theory has also been applied to superconductivity [10], lattice gauge theory [11], and systems of trapped bosons such as Bose–Einstein condensates [12]. These examples and the cited works are by no means exhaustive, but serve to illustrate the flexibility of the CC formalism.

In the early 1980s, Arponen introduced a novel concept into CC theory, namely the bivariational principle [1,13], resulting in the ECC method [1,2,14], and an interpretation of standard CC theory and ECC theory as variational methods in a more general sense, i.e. they are bivariational. However, the ECC method has seen little use in chemistry due to its immense complexity, even for truncated versions. In physics, on the other hand, the ECC model has advantages over standard CC theory that can make it very useful. To illustrate, the ECC method correctly describes symmetry breaking in the Lipkin–Meshkov–Glick quasispin model of collective monopole vibrations in nuclei [13,15], in contrast to the standard CC method, which cannot. For the electronic-structure problem in quantum chemistry, the standard CC model fails dramatically to reproduce dissociation curves of even simple dimers like ${\mathrm{N}}_2$, while the ECC method performs quite well [16,17]. Thus, we are of the opinion that the ECC method is still worthwhile to study, and approximate forms may still prove to be useful in quantum chemistry.

The non-canonical and canonical hierarchies (N)C-ECC(n) introduced in this article turn out to be equivalent, and give identical predictions, when truncated with an excitation-rank complete scheme. On the other hand, the working equations are different and in fact cheaper in the canonical case, albeit marginally. An example is the NC-ECC(1)SD method, i.e. the standard CCSD approach, and the C-ECC(1)SD method, which are equivalent. We also raise the question about diagnostics for practical calculations, and show some numerical evidence that diagnostics can favourably be done using canonical coordinates, even if the computations are done in the usual manner using non-canonical variables.

Another well-known special case is NC-ECC(2), the quadratic coupled-cluster (QCC) method introduced by Van Vorhiis and Head-Gordon [18,19], whose canonical and non-canonical versions are equivalent in their doubles-only approximation. We note that the asymptotic cost of QCCD is the same as CCSD, even if it is a higher-order approximation to ECC. Furthermore, the perfect-pairing (PP) hierarchy [20] of amplitude truncation schemes can be applied to our methods. The PP hierarchy are approximations to the complete-active space self-consistent field (CASSCF) method, including only a tiny subset of even-rank amplitudes combined with orbital optimisation, the latter which we disregard here. The corresponding canonical and non-canonical formulations (N)C-ECC(1)PPH are inequivalent. The n>1 versions could also be interesting in their own right, as investigated by Byrd and coworkers in the case of the QCC method [19].

The remainder of the article is organised as follows: In Section 2, we introduce the bivariational principle and the mathematical setting of local analysis of CC methods. The key concept of our analysis is the notion of local strong monotonicity of the flipped gradient of a smooth bivariational energy functional (see Equation (7(7) $F(T,\mathrm{\Lambda })=({\mathrm{\partial }}_{\mathrm{\Lambda }}\mathcal{E}(T,\mathrm{\Lambda }),{\mathrm{\partial }}_T\mathcal{E}(T,\mathrm{\Lambda })),$(7) )). The usefulness of this property is presented in Theorem 2.1, where local uniqueness and quadratic error estimates are established in a very general setting using Zarantonello's Theorem from nonlinear monotone operator theory [21,22]. Next, Theorem 2.2 summarises the main results of Ref. [23], where strong monotonicity is proven for the non-canonical ECC method. For a recent review on monotonicity in CC theory we refer to [24], where this property is linked to spectral gaps of the systems under study. Section 3 presents the idea of monotonicity-preserving coordinate transformations. Our main result, Theorem 3.1, is a change-of-coordinates result. When combined with Theorems 2.1 and 2.2, the analysis of (N)C-ECC(n) follows in Corollary 3.2. Our tools rely heavily on the functional analytic formulation of cluster operators and the Schrödinger equation developed by Rohwedder and Schneider [25–27]. The results are formulated in a qualitative manner, in the sense that they are indeed rigorous but depend on constants whose numerical (or optimal) values are unknown. We leave further quantitative investigations for future work. In Section 4, we perform some numerical experiments to elucidate some aspects of the (N)C-ECC(n) hierarchies, before we finish with some concluding remarks in Section 5. All the proofs of our results are presented in Appendix.

2. The non-canonical extended coupled-cluster model

2.1. Bivariational principle

The starting point is a generalisation of the Rayleigh–Ritz variational principle to operators that are not necessarily self-adjoint (Hermitian in the finite-dimensional case). For simplicity, we assume a real Hilbert space $\mathcal{H}$. Given a system Hamiltonian $\hat{H}:D(\hat{H})\to \mathcal{H}$, where $D(\hat{H})\subset \mathcal{H}$ is dense, we define a bivariate Rayleigh quotient, ${\mathcal{E}}_{\mathrm{bivar}}:\mathcal{H}\oplus \mathcal{H}\to \mathbb{R}$, (2) ${\mathcal{E}}_{\mathrm{bivar}}(\psi ,\tilde{\psi })=\frac{⟨\tilde{\psi },\hat{H}\psi ⟩}{⟨\tilde{\psi },\psi ⟩},⟨\tilde{\psi },\psi ⟩\ne 0.$(2) Requiring the functional ${\mathcal{E}}_{\mathrm{bivar}}$ to be stationary at $({\psi }_{\ast },{\tilde{\psi }}_{\ast })$ with respect to arbitrary variations in the two wavefunctions leads to the conditions $⟨{\tilde{\psi }}_{\ast },{\psi }_{\ast }⟩\ne 0$, $\hat{H}{\psi }_{\ast }=E_{\ast }{\psi }_{\ast }$ and ${\hat{H}}^{†}{\tilde{\psi }}_{\ast }=E_{\ast }{\tilde{\psi }}_{\ast }$, with $E_{\ast }={\mathcal{E}}_{\mathrm{bivar}}({\psi }_{\ast },{\tilde{\psi }}_{\ast })$, i.e. the right and left eigenvalue problem for $\hat{H}$. If $\hat{H}$ is self-adjoint, the eigenfunctions are identical up to normalisation. The introduction of two independent wavefunctions therefore might seem to complicate matters. However, the bivariate Rayleigh quotient $\mathcal{E}$ allows distinct approximations of ψ and $\tilde{\psi }$, introducing more flexibility for approximate schemes. Moreover, the state defined is a (non-Hermitian) density operator, which is unique, $\rho =\frac{|\psi ⟩⟨\tilde{\psi }|}{⟨\tilde{\psi }|\psi ⟩}.$ When determined variationally, the Hellmann–Feynman theorem [28] gives well-defined physical predictions in terms of ρ.

As is common in analysis of partial differential equations [29,30], we pass to a weak formulation, which in this case is equivalent to the strong formulation outlined above. Under the assumption that $\hat{H}$ is below bounded, we can introduce a unique extension $H:\mathcal{X}\to {\mathcal{X}}^{\prime }$ (dual space), where $\mathcal{X}\subset \mathcal{H}$ is a dense subspace, a Hilbert space with norm $\Vert \cdot {\Vert }_{\mathcal{X}}$, continuously embedded in $\mathcal{H}$. It follows that $\mathcal{H}$ is continuously embedded in ${\mathcal{X}}^{\prime }$, and we have a scale of spaces with dense embeddings, $\mathcal{X}\hookrightarrow \mathcal{H}\hookrightarrow {\mathcal{X}}^{\prime }$. The operator H is bounded (i.e. continuous), and satisfies a Gårding estimate, i.e. for some $\alpha \ge 0$ and some $\mu \in \mathbb{R}$, $⟨\psi ,\mathrm{H\psi }⟩\ge \alpha \Vert \psi {\Vert }_{\mathcal{X}}^2+\mu \Vert \psi {\Vert }^2$ for all $\psi \in \mathcal{X}$. For the electronic-structure problem $\mathcal{X}$ can be taken to be the space of square-integrable functions with finite kinetic energy.

If $\mathcal{H}$ is finite-dimensional, we can set $\mathcal{X}\equiv \mathcal{H}$, simplifying matters a lot, and the reader may if she or he wishes stick to this picture for simplicity, where all operators are basically matrices. In the infinite-dimensional case, however, $\hat{H}$ is typically unbounded as an operator over $\mathcal{H}$, and the above construction is necessary.

Under the stated conditions, ${\mathcal{E}}_{\mathrm{bivar}}:\mathcal{X}\oplus \mathcal{X}\to \mathbb{R}$ is a (Fréchet) smooth map away from the singularity $⟨\tilde{\psi },\psi ⟩=0$, and Taylor series exist and converge locally, allowing a certain degree of intuition to be borrowed from the finite-dimensional case. The right and left Schrödinger equations are then ${\mathrm{\partial }}_{\tilde{\psi }}{\mathcal{E}}_{\mathrm{bivar}}({\psi }_{\ast },{\tilde{\psi }}_{\ast })=0$ and ${\mathrm{\partial }}_{\psi }{\mathcal{E}}_{\mathrm{bivar}}({\psi }_{\ast },{\tilde{\psi }}_{\ast })=0$, respectively. This is the bivariational principle.

2.2. Exponential ansatz and the ECC method

The standard CC method is formulated relative to a fixed reference ${\varphi }_0\in \mathcal{X}$ on determinantal form. By introducing a cluster operator $T=T_1+T_2+\cdots$ with $T_k$ containing all excitations of rank k, i.e. of k fermions relative to ${\varphi }_0$, we have the exact parameterisation $\psi =e^T{\varphi }_0,$ assuming intermediate normalisation, $⟨{\varphi }_0,\psi ⟩=1$.

Since all excitations commute, the cluster operators form a commutative Banach algebra under suitable conditions which we now describe [25,26]. We expand the cluster operators using amplitudes and basis operators, i.e. $T=\sum _{\mu \in \mathcal{I}}{\tau }_{\mu }{X}_{\mu }$, where $X_{\mu }$ excites a number $n=n(\mu )$ of fermions in the reference into the virtual space, i.e. $X_{\mu }=c_{a_1}^{†}{c}_{i_1}\cdots c_{a_n}^{†}{c}_{i_n},$ where the $i_k$ are among the occupied orbitals of ${\varphi }_0$, and $a_k$ among the unoccupied orbitals. The set $\mathcal{I}$ is the generic set of amplitude indices. We introduce a Hilbert space $\mathcal{V}$ with norm $\Vert T\Vert =\Vert T{\varphi }_0{\Vert }_{\mathcal{X}}$, which becomes a useful space for formulating abstract CC theory. Fundamental results include that any $T\in \mathcal{V}$ is a bounded operator on $\mathcal{X}$, such that, e.g. $\exp (T)$ also is a bounded operator. Moreover, $T^{†}$ is also a bounded operator, which means that we can make sense of, e.g. $\exp (-T)\mathrm{Hexp}(T)$, and that we can represent any intermediately normalised $\psi \in \mathcal{X}$ as $\psi =e^T{\varphi }_0$ with $T\in \mathcal{V}$ unique. Finally, all the elements of the algebra are nilpotent. The Banach algebra structure on $\mathcal{V}$ allows CC theory to be rigorously formulated in the full, infinite-dimensional case. This was the approach taken in Ref. [23] for a first analysis of NC-ECC theory.

Again, the finite-dimensional case may by kept in mind: In this case, cluster amplitudes are simply finite-dimensional vectors, and the existence of the exponential parameterisation is a trivial result. There is no need to introduce the norm $\Vert T\Vert$, instead the Euclidean norm on the amplitudes may be used.

Any $\tilde{\psi }$ normalised according to $⟨\tilde{\psi },\psi ⟩=1$ can be represented by introducing a second cluster operator $\mathrm{\Lambda }={\mathrm{\Lambda }}_1+{\mathrm{\Lambda }}_2+\cdots$, viz., $\tilde{\psi }=e^{-T^{†}}{e}^{\mathrm{\Lambda }}{\varphi }_0.$ Inserting the parametrization of ψ and $\tilde{\psi }$ into the bivariate Rayleigh quotient, we obtain the energy functional ${\mathcal{E}}_{\mathrm{ECC}}:\mathcal{V}\oplus \mathcal{V}\to \mathbb{R}$ of the non-canonical ECC method, given in Equation (1(1) ${\mathcal{E}}_{\mathrm{ECC}}(T,\mathrm{\Lambda })=⟨{\varphi }_0,e^{{\mathrm{\Lambda }}^{†}}{e}^{-T}H{e}^T{\varphi }_0⟩,$(1) ). This map is everywhere smooth, and its critical points $(T_{\ast },{\mathrm{\Lambda }}_{\ast })$ are equivalent to the Schrödinger equation and its dual: Under the assumption that the eigenfunctions can be normalised according to $⟨{\varphi }_0,{\psi }_{\ast }⟩=⟨{\tilde{\psi }}_{\ast },{\psi }_{\ast }⟩=1$, ${\psi }_{\ast }$ and ${\tilde{\psi }}_{\ast }$ solve the Schrödinger equation and its dual if and only if (3) $\frac{\mathrm{\partial }{\mathcal{E}}_{\mathrm{ECC}}(T_{\ast },{\mathrm{\Lambda }}_{\ast })}{\mathrm{\partial \Lambda }}=0\mathrm{and}\frac{\mathrm{\partial }{\mathcal{E}}_{\mathrm{ECC}}(T_{\ast },{\mathrm{\Lambda }}_{\ast })}{\mathrm{\partial }T}=0.$(3) Assuming that the eigenvalue $E_{\ast }={\mathcal{E}}_{\mathrm{ECC}}(T_{\ast },{\mathrm{\Lambda }}_{\ast })$ is non-degenerate, $(T_{\ast },{\mathrm{\Lambda }}_{\ast })$ is easily seen to be locally unique.

2.3. Truncations and monotonicity analysis

The non-canonical ECC energy is just one out of many possible parameterisations of the exact bivariate Rayleigh quotient ${\mathcal{E}}_{\mathrm{bivar}}$. In this section, we take a more abstract approach and consider a general energy functional $\mathcal{E}:\mathcal{V}\oplus \mathcal{V}\to \mathbb{R}$, obtained by some exact parameterisation of $(\psi ,\tilde{\psi })$ by means of the space $\mathcal{V}\oplus \mathcal{V}$, i.e. by a pair of cluster operators $(T,\mathrm{\Lambda })$. We will discuss several such functionals in Section 3, obtained from the NC-ECC functional by coordinate transformations.

Only in rare cases can the amplitude equations (3(3) $\frac{\mathrm{\partial }{\mathcal{E}}_{\mathrm{ECC}}(T_{\ast },{\mathrm{\Lambda }}_{\ast })}{\mathrm{\partial \Lambda }}=0\mathrm{and}\frac{\mathrm{\partial }{\mathcal{E}}_{\mathrm{ECC}}(T_{\ast },{\mathrm{\Lambda }}_{\ast })}{\mathrm{\partial }T}=0.$(3) ) be solved exactly. Introduce therefore a discretised space ${\mathcal{V}}_d\subset \mathcal{V}$ of finite dimension by truncating the amplitude index set ${\mathcal{I}}_d\subset \mathcal{I}$, that is, $T_d\in {\mathcal{V}}_d$ if and only if (4) $T_d=\sum _{\mu \in {\mathcal{I}}_d}{\tau }_{d,\mu }{X}_{\mu }\in {\mathcal{V}}_d.$(4) The set ${\mathcal{I}}_d$ is typically defined by the restriction of the excitations to a finite virtual space (a finite basis), and to a finite excitation rank (smaller than the number of electrons). In the chemistry literature, the excitation hierarchy for a given basis is traditionally denoted singles (S), doubles (D), and so on. In the ECC literature, one typically speaks of the SUBn approximation, with n being the maximum rank.

When the discrete space is established, we define a discrete solution by the stationary conditions of the restricted energy function ${\mathcal{E}}_d=\mathcal{E}{\upharpoonright }_{{\mathcal{V}}_d\oplus {\mathcal{V}}_d}$. The stationary equations take the form (5) $\frac{\mathrm{\partial }\mathcal{E}(T_{d\ast },{\mathrm{\Lambda }}_{d\ast })}{\mathrm{\partial }{\lambda }_{\mu }}=\frac{\mathrm{\partial }\mathcal{E}(T_{d\ast },{\mathrm{\Lambda }}_{d\ast })}{\mathrm{\partial }{\tau }_{\mu }}=0,$(5) for all $\mu \in {\mathcal{I}}_d$.

It is not necessary to use the traditional truncation scheme outlined here; any increasing sequence of subspaces ${\mathcal{V}}_d\subset \mathcal{V}$, with d a parameter, that can approximate elements in $\mathcal{V}$ arbitrarily well by increasing d can be used. We let $\mathrm{dist}(v,{\mathcal{V}}_d)$ be the distance from v to ${\mathcal{V}}_d$ measured with respect to the norm of $\mathcal{V}$. Consequently, for all $v\in \mathcal{V}$ we have $\mathrm{dist}(v,{\mathcal{V}}_d)\to 0$ as $d\to +\mathrm{\infty }$. Such a sequence of spaces is referred to as a Galerkin sequence. Other options than the traditional truncation schemes are explicitly correlated methods [31] and complete-active space methods [20,32,33] such as the PP hierarchy.

An often overlooked point in the physics literature is the fact that convergence of the equations does not in general imply convergence of their solutions. An important question is therefore whether the discrete critical points $(T_{d\ast },{\mathrm{\Lambda }}_{d\ast })$ converge to the exact critical points $(T_{\ast },{\mathrm{\Lambda }}_{\ast })$ as $d\to +\mathrm{\infty }$. This would imply that the energy converges too, and in a quadratic manner due to the critical point formulation.

Monotonicity is an important notion in connection with the local analysis of the CC method and its variations [23–27,34]. The use of monotonicity in the analysis of the standard CC method was introduced by Schneider and Rohwedder [25,27]. The intuition behind (strict) monotonicity is that of an everywhere (strictly) increasing (or decreasing) function. The function that one studies is typically a root problem such as the CC amplitude equations $f(T)=0$, where $f_{\mu }(T)=⟨{\varphi }_{\mu },e^{-T}H{e}^T{\varphi }_0⟩$. Monotonicity allows the establishment of locally unique solutions of the Galerkin problem and is therefore important for the motivation of numerical implementations. As such it is a fundamental result of the CC method's practical usage in quantum chemistry. It also connects spectral gaps, e.g. HOMO-LUMO gap, to stability constants within the analysis [24]. (See also the steerable CAS-ext gap connected to the tailored CC method [35] that treats quasi-degenerate systems [36].)

The particular monotonicity property that is key for this presentation is an even stronger version than that of strict monotonicity, and is called strong monotonicity. It is defined as follows: A finite-dimensional vector-valued function $F(Z)$ is locally strongly monotone near some $Z_{\ast }$ if for $Z_1,Z_2$ in a neighbourhood of $Z_{\ast }$ we have (6) $⟨F(Z_1)-F(Z_2),Z_1-Z_2⟩\ge \eta \Vert Z_1-Z_2{\Vert }^2,$(6) for some constant $\eta >0$. (In the infinite-dimensional case $⟨\cdot ,\cdot ⟩$ is the dual pairing, which then becomes an infinite sum, see Ref. [23] for more details.)

Furthermore, we need the concept of Lipschitz continuity: F is locally Lipschitz continuous with constant L>0 if $\Vert F(Z_1)-F(Z_2)\Vert \le L\Vert Z_1-Z_2\Vert .$ In particular, any (Fréchet) smooth function is locally Lipschitz continuous, and so are all its derivatives.

The map F that we will study is the flipped gradient of the general energy functional $\mathcal{E}:\mathcal{V}\oplus \mathcal{V}\to \mathbb{R}$, defined as (7) $F(T,\mathrm{\Lambda })=({\mathrm{\partial }}_{\mathrm{\Lambda }}\mathcal{E}(T,\mathrm{\Lambda }),{\mathrm{\partial }}_T\mathcal{E}(T,\mathrm{\Lambda })),$(7) or more compactly $F(T,\mathrm{\Lambda })=R\mathrm{\partial }\mathcal{E}(T,\mathrm{\Lambda })$, with R being the map that exchanges the partial derivatives. The motivation is as follows: If we consider the bivariate Rayleigh quotient, $\mathrm{\partial }{\mathcal{E}}_{\mathrm{bivar}}$ is not locally strongly monotone, as its critical points are saddle points. On the other hand, the flipped gradient $F_{\mathrm{bivar}}=R\mathrm{\partial }{\mathcal{E}}_{\mathrm{bivar}}$ can be seen to be locally strongly monotone near the ground state, given that this ground state is non-degenerate with a nonzero spectral gap to the remaining spectrum. It is natural to expect that one can find conditions such that the flipped gradient of the energy when expressed in new coordinates is locally strongly monotone.

The following is a central result, adapting a result due to Zarantonello [21,22] (points 1 and 2) to the present notation and setting, and applied to the flipped gradient of an energy functional (point 3).

Theorem 2.1

General convergence and error estimates

Let $F:\mathcal{V}\oplus \mathcal{V}\to {\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime }$ be a map, and let $U\subset \mathcal{V}\oplus \mathcal{V}$ be an open ball containing a $Z_{\ast }$ such that $F(Z_{\ast })=0$.

Let ${\mathcal{V}}_d\subset \mathcal{V}$ be a Galerkin sequence of subspaces with $P_d$ being the orthogonal projector onto ${\mathcal{V}}_d\oplus {\mathcal{V}}_d$. Furthermore, let $F_d:{\mathcal{V}}_d\oplus {\mathcal{V}}_d\to {\mathcal{V}}_d^{\prime }\oplus {\mathcal{V}}_d^{\prime }$ be the Galerkin discretisation of F, i.e. $F_d(Z_d)=P_{d}F(Z_d)$.

Assume that F is locally strongly monotone with constant $\eta >0$ and Lipschitz continuous with constant L>0 on U. Then, the following holds:

(1)	$Z_{\ast }$ is the only root in U.
(2)	There is a sufficiently large $d_0$, such that for any $d>d_0$, there exists $Z_{d\ast }\in {\mathcal{V}}_d\oplus {\mathcal{V}}_d$ such that $F_d(Z_{d\ast })=0$. This root is unique in U and we have the following error estimate (quasi-optimality of the discrete solution): (8) $\Vert Z_{\ast d}-Z_{\ast }\Vert \le \frac{L}{\eta }\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d).$(8)

Let $\mathcal{E}:\mathcal{V}\oplus \mathcal{V}\to \mathbb{R}$, $Z\mapsto \mathcal{E}(Z)$ be a (Fréchet) smooth energy functional. Let R be the flipping map as introduced after Equation (7(7) $F(T,\mathrm{\Lambda })=({\mathrm{\partial }}_{\mathrm{\Lambda }}\mathcal{E}(T,\mathrm{\Lambda }),{\mathrm{\partial }}_T\mathcal{E}(T,\mathrm{\Lambda })),$(7) ) and set $F=R\mathrm{\partial }\mathcal{E}$, and $E_{\ast }=\mathcal{E}(Z_{\ast })$.

(3)

For $d>d_0$, the discrete Galerkin equations $\mathrm{\partial }{\mathcal{E}}_d(Z_{\ast d})=0$ have locally unique solutions, and in addition to the error estimate (8(8) $\Vert Z_{\ast d}-Z_{\ast }\Vert \le \frac{L}{\eta }\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d).$(8) ), we have the energy error (9) $\begin{array}{rl}|\mathcal{E}(Z_{\ast d})-E_{\ast }|& \le C\Vert Z_{\ast d}-Z_{\ast }{\Vert }^2\\ & \le C{\left(\frac{L}{\eta }\right)}^2\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d{)}^2.\end{array}$(9)

The proof is presented in Appendix. The error estimate (9(9) $\begin{array}{rl}|\mathcal{E}(Z_{\ast d})-E_{\ast }|& \le C\Vert Z_{\ast d}-Z_{\ast }{\Vert }^2\\ & \le C{\left(\frac{L}{\eta }\right)}^2\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d{)}^2.\end{array}$(9) ) shows that for (smooth) energy functionals with a locally strongly monotone flipped gradient, the bivariational method of discretisation behaves very similar to the usual Rayleigh–Ritz variational method of discretisation. As we enlarge the Galerkin space, the discrete ground state converges, and the energy error is quadratic in the error of the state. However, we cannot guarantee convergence from above, but this is much less important than actually having a quadratic error.

An interesting fact is that Brouwer's fixed point theorem [37] can be used to obtain a sufficient condition for the constant $d_0$, where quadratic convergence sets in, namely (10) $\mathrm{dist}(Z_{\ast },{\mathcal{V}}_{d_0}\oplus {\mathcal{V}}_{d_0})<\frac{\delta }{1+L/\eta },$(10) see Refs. [23,26]. The radius δ of the domain U in Theorem 2.1 is unknown in general, and we see that a small monotonicity constant η relative to L will be severely detrimental for the convergence, as it forces $d_0$ to be very large. On the other hand, in general one has that η and the Lipschitz constant L are related by $\eta \le L$. The optimal value of the right-hand side is therefore $\delta /2$.

The following theorem summarises the main results of Ref. [23], where the proof and more details can be found:

Theorem 2.2

NC-ECC monotonicity

Assume that the system Hamiltonian $\hat{H}$ is self-adjoint, and that the ground state of $\hat{H}$ exists, is non-degenerate, and that there is a spectral gap $\gamma >0$ between the ground-state energy $E_{\ast }$ and the rest of the spectrum. Assume that the reference ${\varphi }_0$ is such that it is not orthogonal to the ground-state wavefunction. Let $Z_{\ast }=(T_{\ast },{\mathrm{\Lambda }}_{\ast })\in \mathcal{V}\oplus \mathcal{V}$ be the corresponding critical point of ${\mathcal{E}}_{\mathrm{ECC}}$, and assume that $T_{\ast }$ and ${\mathrm{\Lambda }}_{\ast }$ are not too large, i.e. that ${\varphi }_0$ is a sufficiently good approximation to ${\psi }_{\ast }$. Then, $F=R\mathrm{\partial }{\mathcal{E}}_{\mathrm{ECC}}$ is locally strongly monotone near $Z_{\ast }$ with a constant $\eta =\mathrm{C\gamma }$, for some C<1.

The consequence of Theorem 2.2, when combined with Theorem 2.1, is that the non-canonical ECC method is convergent as the discrete cluster amplitude space ${\mathcal{V}}_d$ approaches the untruncated limit. As already remarked, we do not explicitly know the onset $d_0$ of quadratic convergence.

3. Monotonicity-preserving coordinate transformations

3.1. A class of exact coupled-cluster models

In addition to the non-canonical ECC parameterisation, Arponen also considered a second parameterisation of the bra and ket wavefunctions, which gives equations of motion for the time-dependent Schrödinger equation that are canonical in the sense of Hamiltonian mechanics [1,2]. (This must not be confused with the use of canonical Hartree–Fock orbitals, which is unrelated.) This parameterisation is given in terms of a coordinate transformation ${\theta }_{\mathrm{C}-\mathrm{ECC}}:\mathcal{V}\oplus \mathcal{V}\to \mathcal{V}\oplus \mathcal{V}$ as (11) $(T,\mathrm{\Lambda })={\theta }_{\mathrm{C}-\mathrm{ECC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime }),$(11) where ${\mathrm{\Lambda }}^{\prime }=\mathrm{\Lambda }$, and $T=S(T^{\prime };{\mathrm{\Lambda }}^{\prime })$ is defined by (12) $\mathrm{QT}{\varphi }_0=Q{e}^{-{{\mathrm{\Lambda }}^{\prime }}^{†}}{T}^{\prime }{\varphi }_0,Q=I-|{\varphi }_0⟩⟨{\varphi }_0|.$(12) This function has inverse $Q{T}^{\prime }{\varphi }_0=Q{e}^{{\tilde{\mathrm{\Lambda }}}^{†}}T{\varphi }_0$. (In Arponen's work [2], the notation $(T^{\prime },{\mathrm{\Lambda }}^{\prime })=(\mathrm{\Sigma },{\tilde{\mathrm{\Sigma }}}^{†})$ is used.) The map ${\theta }_{\mathrm{C}-\mathrm{ECC}}$ is smooth and invertible with a smooth inverse, and we therefore obtain a new exact energy functional ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}}={\mathcal{E}}_{\mathrm{ECC}}\circ {\theta }_{\mathrm{C}-\mathrm{ECC}},$ with values (13) ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=⟨{\varphi }_0,e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{e}^{-S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}H{e}^{S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}{\varphi }_0⟩.$(13) A remarkable consequence of this second parameterisation is that it corresponds to retaining only those terms in Equation (1(1) ${\mathcal{E}}_{\mathrm{ECC}}(T,\mathrm{\Lambda })=⟨{\varphi }_0,e^{{\mathrm{\Lambda }}^{†}}{e}^{-T}H{e}^T{\varphi }_0⟩,$(1) ) that can be represented by ‘doubly linked’ diagrams [1,2], (14) ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })={⟨{\varphi }_0,e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}(H{e}^{T^{\prime }}{)}_{\mathrm{C}}{\varphi }_0⟩}_{\mathrm{DL}}.$(14) The phrase ‘doubly linked’ means that every power of $({\mathrm{\Lambda }}^{\prime }{)}^{†}$ is connected to two $T^{\prime }$ operators on its right, unless it is connected directly to H. The subscript ‘C’ for ‘connected’ is the usual connectedness criterion on contractions between H and powers of $T^{\prime }$ [8]. Thus, the canonical coordinates represent a more compact representation in that the resulting tensor contractions or diagrams in the energy are identical to those obtained in the NC-ECC energy (1(1) ${\mathcal{E}}_{\mathrm{ECC}}(T,\mathrm{\Lambda })=⟨{\varphi }_0,e^{{\mathrm{\Lambda }}^{†}}{e}^{-T}H{e}^T{\varphi }_0⟩,$(1) ), except for some diagrams that are explicitly eliminated.

Similarly, for the standard CC method, Arponen introduced the coordinate transformation ${\theta }_{\mathrm{CC}}$ given by (15) $(T,\mathrm{\Lambda })={\theta }_{\mathrm{CC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=(T^{\prime },e^{{\mathrm{\Lambda }}^{\prime }}-1).$(15) We obtain the energy functional ${\mathcal{E}}_{\mathrm{CC}}={\mathcal{E}}_{\mathrm{ECC}}\circ {\theta }_{\mathrm{CC}}$, where (16) ${\mathcal{E}}_{\mathrm{CC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=⟨{\varphi }_0,(1+({\mathrm{\Lambda }}^{\prime }{)}^{†})e^{-T^{\prime }}H{e}^{T^{\prime }}{\varphi }_0⟩,$(16) that is, the standard CC Lagrangian [3]. Incidentally, the standard CC coordinates are also canonical.

The map ${\theta }_{\mathrm{CC}}$ can be generalised to Taylor polynomials. By setting $e^{\mathrm{\Lambda }}=(e^{{\mathrm{\Lambda }}^{\prime }}{)}_n\equiv 1+{\mathrm{\Lambda }}^{\prime }+\frac{1}{2}({\mathrm{\Lambda }}^{\prime }{)}^2+\cdots +\frac{1}{n!}({\mathrm{\Lambda }}^{\prime }{)}^n,$ we can solve for Λ in terms of ${\mathrm{\Lambda }}^{\prime }$ by, e.g. considering first the singles, then doubles, etc., giving a smooth map $G_n:\mathcal{V}\to \mathcal{V}$ such that $e^{G_n({\mathrm{\Lambda }}^{\prime })}=(e^{{\mathrm{\Lambda }}^{\prime }}{)}_n$. In fact, since the cluster operators are nilpotent, $G_n({\mathrm{\Lambda }}^{\prime })=\ln [(e^{{\mathrm{\Lambda }}^{\prime }}{)}_n]$, where the logarithm is expanded in a (finite) Taylor series around the identity. Similarly, we can solve for ${\mathrm{\Lambda }}^{\prime }$ in terms of Λ, demonstrating that this map has an inverse, and in fact that this inverse is smooth. We obtain a coordinate transformation ${\theta }_n$ given by (17) $(T,\mathrm{\Lambda })={\theta }_n(T^{\prime },{\mathrm{\Lambda }}^{\prime })=(T^{\prime },G_n({\mathrm{\Lambda }}^{\prime })),$(17) and the corresponding energy functional (18a) ${\mathcal{E}}_{\mathrm{NC}-\mathrm{ECC}(n)}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n{e}^{-T^{\prime }}H{e}^{T^{\prime }}{\varphi }_0⟩.$(18a) Coordinate transformations form a group, and may thus be composed. By combining ${\theta }_{\mathrm{C}-\mathrm{ECC}(n)}={\theta }_n\circ {\theta }_{\mathrm{C}-\mathrm{ECC}}$, we obtain an energy functional (18b) $\begin{array}{rl}& {\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}(n)}(T^{\prime },{\mathrm{\Lambda }}^{\prime })\\ & =⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n{e}^{-S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}H{e}^{S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}{\varphi }_0⟩.\end{array}$(18b) Since, in the NC-ECC energy functional (1(1) ${\mathcal{E}}_{\mathrm{ECC}}(T,\mathrm{\Lambda })=⟨{\varphi }_0,e^{{\mathrm{\Lambda }}^{†}}{e}^{-T}H{e}^T{\varphi }_0⟩,$(1) ), an exponential $e^{-{\mathrm{\Lambda }}^{†}}$ can be inserted after $e^T$ without changing the result, both of these hierarchies correspond to truncations of a Baker–Campbell–Hausdorff expansion at order n, and are thus manifestly extensive.

3.2. Coordinate transformation theorem

Equations (18a(18a) ${\mathcal{E}}_{\mathrm{NC}-\mathrm{ECC}(n)}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n{e}^{-T^{\prime }}H{e}^{T^{\prime }}{\varphi }_0⟩.$(18a) ) and (18b(18b) $\begin{array}{rl}& {\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}(n)}(T^{\prime },{\mathrm{\Lambda }}^{\prime })\\ & =⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n{e}^{-S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}H{e}^{S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}{\varphi }_0⟩.\end{array}$(18b) ) represent two hierarchies of exact parameterisations of the bivariate Rayleigh quotient. It is therefore of interest to determine whether they have locally strongly monotone flipped gradients. To establish this, we study the effect on local strong monotonicty of a coordinate transformation.

Theorem 3.1

Coordinate transformations

Let $\mathcal{E}:\mathcal{V}\oplus \mathcal{V}\to \mathbb{R}$ be a smooth energy functional, let $Z_{\ast }$ be a critical point, and assume that $F=R\mathrm{\partial }\mathcal{E}$ is locally strongly monotone near $Z_{\ast }$ with constant $\eta >0$. Let a smooth $\theta :\mathcal{V}\oplus \mathcal{V}\to \mathcal{V}\oplus \mathcal{V}$ with a smooth inverse be a given coordinate transformation, and let ${\mathcal{E}}_{\theta }=\mathcal{E}\circ \theta$ be the energy functional expressed in the new coordinates. Let $W_{\ast }={\theta }^{-1}(Z_{\ast })$ be the corresponding critical point for ${\mathcal{E}}_{\theta }$, and let $F_{\theta }=R\mathrm{\partial }{\mathcal{E}}_{\theta }$ be its flipped gradient. Let $M_{\ast }=\mathrm{\partial }\theta (W_{\ast })$ be the Jacobian at $W_{\ast }$. Then we have the following conclusions:

(1)	If $M_{\ast }R=R{M}_{\ast }$, then $F_{\theta }$ is locally strongly monotone near $W_{\ast }$ with constant $\Vert M_{\ast }^{-1}{\Vert }^{-2}\eta$.
(2)	In the noncommuting case, if $m_{\ast }=M_{\ast }-I$ is sufficiently small, $F_{\theta }$ is locally strongly monotone near $W_{\ast }$ with constant ${\eta }^{\prime }=\eta \Vert M_{\ast }^{-1}{\Vert }^{-2}-C(I+\Vert m_{\ast }\Vert )\Vert m_{\ast }\Vert ,$where C is the constant from Theorem 2.1(3).

Theorem 3.1 tells us that if we create a new method by changing coordinates in a sufficiently well-behaved manner, the monotonicity of the flipped gradient of the energy will be preserved. Hence, the new method will be convergent in the sense of Theorem 2.1, if the original method is. This gives us great flexibility in choosing new cluster-operator parameterisations, since we already know that non-canonical ECC is convergent. The interested reader can find the proof in the Appendix.

3.3. Monotonicity of (N)C-ECC$(n)$ models

We apply Theorem 3.1 to the maps ${\theta }_n$ and ${\theta }_n\circ {\theta }_{\mathrm{C}-\mathrm{ECC}}$ that define the NC-ECC$(n)$ and C-ECC$(n)$ models, respectively. The conclusion is as follows (with proof given in the Appendix):

Corollary 3.2

For any of the NC-ECC$(n)$ or C-ECC$(n)$ models, the assumption that the ground-state critical point $W_{\ast }=(T_{\ast }^{\prime },{\mathrm{\Lambda }}_{\ast }^{\prime })$ is not too large together with a spectral gap $\gamma >0$ is sufficient to guarantee local strong monotonicity of the flipped gradient of the energy, and hence a quasi-optimal solution to the Galerkin problem and a quadratic error estimate for the energy.

We note that our estimates are probably pessimistic for some of the models covered here. The analysis starts with a given monotonicity constant η for the NC-ECC scheme, depending on the spectral gap, and consistently produces an ${\eta }^{\prime }<\eta$ for the method obtained using the coordinate change, worsening the error estimates. It may well be that a direct analysis of the secondary method yields a better ${\eta }^{\prime }$. The assumptions on the Hamiltonian and reference may also become milder. However, the important point here is that Theorem 3.1 does guarantee that the new method is convergent under some reasonable conditions. For example, we have now proven that QCC theory [18] is convergent if the reference is a sufficiently good approximation to the ground state, and using Equation (10(10) $\mathrm{dist}(Z_{\ast },{\mathcal{V}}_{d_0}\oplus {\mathcal{V}}_{d_0})<\frac{\delta }{1+L/\eta },$(10) ) also a basic means to study which truncations or Galerkin schemes can be reasonable, at least in principle. It may be interesting to see whether truncation schemes like the PP hierarchy with orbital optimisation, also in a quadratic n = 2 or higher formulation [19], can be further analysed based on our results.

In the proof of Theorem 3.1 (see Appendix), it arises naturally that the most favourable coordinate transformations are those that commute with the flipping map R, since local strong monotonicity then follows with no assumptions on the Jacobian of the map. The Jacobian commutes with R if and only if (19) $\frac{\mathrm{\partial }T}{\mathrm{\partial }{T}^{\prime }}=\frac{\mathrm{\partial \Lambda }}{\mathrm{\partial }{\mathrm{\Lambda }}^{\prime }},\mathrm{and}\frac{\mathrm{\partial }T}{\mathrm{\partial }{\mathrm{\Lambda }}^{\prime }}=\frac{\mathrm{\partial \Lambda }}{\mathrm{\partial }{T}^{\prime }},$(19) and one must assume that this holds at every point in $\mathcal{V}\oplus \mathcal{V}$, as one does not know a priori where the critical point is. It is not clear what such transformations in general look like, and whether such transformations are useful reparameterisations of the energy.

3.4. Properties of the canonical and non-canonical schemes

The coordinate transformation ${\theta }_{\mathrm{C}-\mathrm{ECC}}$ as represented by Equation (A7(A7) $A({\mathrm{\Lambda }}^{\prime })T^{\prime }=\sum _{\mu ,\nu \in \mathcal{I}}{X}_{\mu }⟨{\varphi }_{\mu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\nu }⟩⟨{\varphi }_{\nu },T^{\prime }{\varphi }_0⟩.$(A7) ) is such that when applied to a cluster operator $T_k^{\prime }$ of rank k, it generates terms $T_{k^{\prime }}$ with $k^{\prime }\le k$. The same is true for the inverse map. Thus, if ${\mathcal{V}}_d$ is excitation-rank complete, i.e. it contains all excitations of rank up to and including some k>0, then the Galerkin discretisation of the NC-ECC(n) method commutes with changing coordinates via the coordinate map ${\theta }_{\mathrm{C}-\mathrm{ECC}}$, ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}(\mathrm{n}),d}={\mathcal{E}}_{\mathrm{NC}-\mathrm{ECC}(\mathrm{n}),d}\circ {\theta }_{\mathrm{C}-\mathrm{ECC}}.$ By inspection, one can also see that this holds for a doubles-only truncation, since $A({\mathrm{\Lambda }}^{\prime })=I$ in this case. We obtain the following result:

Theorem 3.3

Equivalence of canonical and non-canonical coordinates

Let ${\mathcal{V}}_d$ be excitation-rank complete or consist of doubles excitations only, and let n be given. Then, the discrete solutions $(T_{\ast d},{\mathrm{\Lambda }}_{\ast d})$ and $(T_{\ast d}^{\prime },{\mathrm{\Lambda }}_{\ast d}^{\prime })$ of the NC-ECC(n) and C-ECC(n) methods, respectively, are equivalent and related via ${\theta }_{\mathrm{C}-\mathrm{ECC}}$, i.e. $T_{\ast d}=A({\mathrm{\Lambda }}_{\ast d}^{\prime })T_{\ast d}^{\prime }$ and ${\mathrm{\Lambda }}_{\ast d}={\mathrm{\Lambda }}_{\ast d}^{\prime }$.

We stress that, if ${\mathcal{V}}_d$ is not excitation-rank complete, the canonical and non-canonical parameterisations are not equivalent. This would be the case for the PP hierarchy of truncations [19,20].

According to the doubly linked structure of the energy functional ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}(\mathrm{n})}$, see the discussion after Equation (14(14) ${\mathcal{E}}_{\mathrm{C}-\mathrm{ECC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })={⟨{\varphi }_0,e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}(H{e}^{T^{\prime }}{)}_{\mathrm{C}}{\varphi }_0⟩}_{\mathrm{DL}}.$(14) ), the amplitude equations for the canonical case are cheaper, albeit by a small amount. Moreover, it is reasonable to expect that the canonical solution $(T_{\ast d}^{\prime },{\mathrm{\Lambda }}_{\ast d}^{\prime })$ is more compact compared to $(T_{\ast d},{\mathrm{\Lambda }}_{\ast d})$. We investigate this claim numerically in Section 4.

4. Numerical results

4.1. Implementation

The (un-)truncated (N)C-ECC(n)SD Equations (18(18a) ${\mathcal{E}}_{\mathrm{NC}-\mathrm{ECC}(n)}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n{e}^{-T^{\prime }}H{e}^{T^{\prime }}{\varphi }_0⟩.$(18a) ) and (5(5) $\frac{\mathrm{\partial }\mathcal{E}(T_{d\ast },{\mathrm{\Lambda }}_{d\ast })}{\mathrm{\partial }{\lambda }_{\mu }}=\frac{\mathrm{\partial }\mathcal{E}(T_{d\ast },{\mathrm{\Lambda }}_{d\ast })}{\mathrm{\partial }{\tau }_{\mu }}=0,$(5) ), together with the untruncated coordinate transformation Equation (A7(A7) $A({\mathrm{\Lambda }}^{\prime })T^{\prime }=\sum _{\mu ,\nu \in \mathcal{I}}{X}_{\mu }⟨{\varphi }_{\mu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\nu }⟩⟨{\varphi }_{\nu },T^{\prime }{\varphi }_0⟩.$(A7) ) have been implemented in a local full CI-based program, i.e. all intermediates are expressed as vectors in the full CI basis. To this end, the C-ECC(n) amplitudes are computed using transformed NC-ECC(n) residual expressions [17]. The C-ECC(n) amplitude equations are thus: (20a) $\begin{array}{rl}0& =\sum _{\nu \in \mathcal{I}}⟨{\varphi }_{\nu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\mu }⟩⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n[H^S,X_{\nu }]{\varphi }_0⟩,\end{array}$(20a) (20b) $\begin{array}{rl}0& =⟨{\varphi }_{\mu },(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_{n-1}{H}^S{\varphi }_0⟩\\ & -\sum _{\nu \in \mathcal{I}}⟨{\varphi }_{\mu },X_{\nu }^{†}{e}^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{T}^{\prime }{\varphi }_0⟩⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n[H^S,X_{\nu }]{\varphi }_0⟩,\end{array}$(20b) where $H^S=e^{-S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}H{e}^{S(T^{\prime };{\mathrm{\Lambda }}^{\prime })}$. The sparsity of the coordinate transformation Equation (A7(A7) $A({\mathrm{\Lambda }}^{\prime })T^{\prime }=\sum _{\mu ,\nu \in \mathcal{I}}{X}_{\mu }⟨{\varphi }_{\mu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\nu }⟩⟨{\varphi }_{\nu },T^{\prime }{\varphi }_0⟩.$(A7) ) has been exploited throughout, e.g. only singles amplitudes are transformed in the case where ${\mathcal{V}}_d$ contains only singles and doubles. The coupled amplitude equations (20(20a) $\begin{array}{rl}0& =\sum _{\nu \in \mathcal{I}}⟨{\varphi }_{\nu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\mu }⟩⟨{\varphi }_0,(e^{({\mathrm{\Lambda }}^{\prime }{)}^{†}}{)}_n[H^S,X_{\nu }]{\varphi }_0⟩,\end{array}$(20a) ) are solved iteratively starting from an MP2-guess, using an alternating scheme and applying DIIS convergence acceleration. In all computations, residuals and energies were converged to a threshold of ${10}^{-4}$ and ${10}^{-6}$ a.u., respectively. The (N)C-ECC(1)SD and (N)C-ECC(∞)SD implementations are verified by reproducing the ‘CCSD’ and ‘ECCSD’ energies presented in [17].

4.2. Numerical experiments

The (N)C-ECC(n)SD and (N)C-ECC(n)DT models have been studied numerically by investigating the potential energy curves of the hydrogen fluoride molecule with interatomic distances $1.0\le R_{\mathrm{H}-\mathrm{F}}\le 3.5(a_0)$ in a DZV basis set [17] and the nitrogen molecule with $1.5\le R_{\mathrm{N}-\mathrm{N}}\le 4.0(a_0)$ with a frozen core in a 6-31G basis [38]. Additionally, the H${}_8$ model system with structural parameters $0.0001\le \alpha \le 1.0$ in a minimal basis set [32,39] is considered. For large distances R and small α, respectively, these systems comprise significant multireference character, i.e. the weight of the Hartree–Fock configuration in the full CI wave function, $|⟨{\varphi }_0,{\psi }^{\mathrm{FCI}}⟩{|}^2$, is fairly small. Thus, these species are good candidates to study novel quantum chemical methods.

The energy curves of the canonical models C-ECC(n)SD are identical to the non-canonical NC-ECC(n)SD ones and are thus not presented here. However, the results differ if excitation-rank incomplete truncation schemes are employed. For instance, in canonical ECC(n)DT, singles amplitudes are effectively generated from doubles and triples amplitudes, while these are absent in the non-canonical model. This effect has been studied on the potential curve of the HF molecule and is depicted in Figure 1: The generation of singles amplitudes entails that the canonical computation is lower in energy, in particular towards the multireference region where these contribute significantly to the wave function expansion. This depends, however, on the role of the singles amplitudes in the wave function: In test computations on the H${}_8$ model system a different trend was observed, consistent with the diminished importance of singles in the wave function (vide infra). Therefore, we cannot conclude that the canonical coordinates are consistently better when unconventional truncations are used.

Figure 1. Difference between a canonical and a non-canonical ECC(n)DT computation for the potential curve of HF.

In order to investigate the effect of using different coordinates in ECC(n)SD computations, we calculated a set of CC diagnostics which are often used to assess the quality of CC computations [40,41]. These are based on the largest singular value (${\mathcal{D}}_1$) and Frobenius norm (${\mathcal{T}}_1$) of the matrix representation of the singles amplitudes. (Equivalently, ${\mathcal{T}}_1^2=\Vert {\tau }_{\ast 1}{\Vert }_2^2/N$, the sum of the squares of the singles amplitudes, with N the number of correlated electrons.) Although diagnostics based on doubles amplitudes are preferred, they are not as available in implementations as are the singles-based variants [42]. Additionally, we computed the diagnostic $||{\tau }_{\ast }|{|}_2^2/(||{\tau }_{\ast }|{|}_2^2+1)$ which involves all the amplitudes. This choice can be motivated from monotonicity arguments and will be discussed in a forthcoming paper.

The values have been computed for truncation schemes $n=1,2,\mathrm{\infty }$. Since the trends are very similar, only the data for n = 1 is presented for HF and H${}_8$. For N${}_2$, we focus on the results for n = 2, because the quadratic model is well-known to be a better approximation to full CI compared to standard CCSD in this case [16] (cf. line starting upper left with triangle markers in Figure 4). We stress, that even if the NC-ECC(1)SD and C-ECC(1)SD (which is standard CCSD) are equivalent, the singles amplitudes in their respective parameterisations differ, producing different values for the diagnostics. Figure 2 shows the diagnostics correlated with the multireference character for the HF potential curve, in Figures 3 and 4 values for the H${}_8$ model and the N${}_2$ molecule are shown, respectively. In both H${}_8$ and N${}_2$, electron correlation is dominated by doubles amplitudes, as can be seen from the small values of the singles-based diagnostics. Since the reparameterisation of the amplitudes in the canonical model does not affect the amplitudes of highest excitation rank, the difference between the NC-ECC(n)SD and C-ECC(n)SD amplitude vectors is negligible. This is different for the HF case. Here, the amplitude norms of the canonical models are consistently smaller than the non-canonical variants, indicating that the wave function parameterisation is more compact.

Figure 2. Comparison of CC diagnostics of the C-ECC(1)SD and NC-ECC(1)SD model for the HF potential curve correlated with the multireference character.

Figure 3. Comparison of CC diagnostics of the C-ECC(1)SD and NC-ECC(1)SD model for the H${}_8$ potential curve correlated with the multireference character.

Figure 4. Comparison of CC diagnostics of the C-ECC(2)SD and NC-ECC(2)SD model (lower set of five curves, almost on top of each other) for the N${}_2$ potential curve correlated with the multireference character. The red curve shows diagnostics for NC-ECC(1)SD, indicating the inferiority of this model in the multireference region.

Our numerical experiments suggest, that for excitation-rank incomplete models, the canonical map generates effectively an excitation-rank complete parameterisation, but does not necessarily yield significantly better results. Concerning the a priori excitation-rank complete models, it has been found that the canonical parameterisation can be more compact compared to the non-canonical one, a desired property for post-Hartree–Fock methods. In particular, it may be useful to consider using canonical coordinates for diagnostics, even for well-established numerical codes using standard CC formalism.

5. Concluding remarks

In this article, we formulated basic error estimates for a class of exact models, defined in terms of replacing, in Arponen's ECC method, the exact exponential $e^{{\mathrm{\Lambda }}^{†}}$ of the dual cluster operator with an n-th Taylor polynomial, the canonical C-NCC(n) models and the non-canonical NC-ECC(n) models. The central result was a coordinate-transformation theorem, Theorem 3.1, that gives error estimates for any method that can be described as a coordinate transformation of ECC theory. Notably, these results guarantee asymptotically quadratic error estimates for the ground-state energy of all models, under certain mild conditions.

Apart from Theorem 3.1, a basically self-contained mathematical framework for local error analysis of coupled-cluster methods was presented. This was based on Arponen's bivariational principle and basic results from nonlinear monotone operator theory, i.e. Zarantonello's theorem. Also central was our prior analysis of Arponen's extended coupled-cluster method in its non-canonical formulation.

The methods covered by our analysis include standard CC theory, quadratic CC theory [18,19], the perfect-pairing hierarchy [20] for approximating CASSCF, also in its quadratic version [19], and, Arponen's canonical ECC method. The computational cost of the (N)C-ECC$(n)$ methods truncated with the standard singles, doubles, etc., scheme, are expensive: already the nonlinear terms in $({\mathrm{\Lambda }}^{†}{)}^2$ pushes the cost beyond standard CC theory when higher than doubles are considered. However, having obtained exact mathematical characterisations of the hierarchies of methods is an important first step in producing cheaper and more reliable methods compared to standard CC. One such approach could be similar to the CC2 method of Christiansen and coworkers [43], where the cost of CCSD is reduced from $\mathcal{O}(N^6)$ to $\mathcal{O}(N^5)$, $N$ being the system size.

The error estimates are not optimal for many methods. A direct analysis of canonical ECC would probably provide the most optimistic analysis for all the C-ECC(n) methods, due to the doubly linked structure and the equivalence of excitation-rank complete Galerkin discretisations.

Finally, we performed some simple numerical experiments, focusing on the possibility of using canonical coordinates in place of the usual CC amplitudes when doing diagnostic estimates on CC calculations on systems with multireference character. Our preliminary findings support the hypothesis that the canonical coordinates are more compact compared to the usual coordinates, providing more accurate diagnostics.

An interesting extension of the present work would be to study truncations where singles-amplitudes are replaced by orbital rotations, either unitary or biorthogonal, as in the QCC and PP approaches, or the non-orthogonal orbital optimised coupled-cluster method of Pedersen and coworkers. [44]. Moreover, the complete-active space coupled-cluster method by Adamowicz and coworkers [32] fits the present scheme. It is also known that quadratic CC and ECC in general are quite good at reproducing multireference character, while standard single-reference CC is quite poor at this. Thus, a modified analysis of the ECC method that includes multireference assumptions, such as the steerable CAS-ext gap of Ref. [34], could potentially lead to a deeper understanding of how CC methods generally behave in the presence of static correlation.

Acknowledgments

The authors are thankful to M. A. Csirik for useful comments.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This work has received funding from the Research Council of Norway (RCN) [CoE grant numbers 287906 and 262695] (Hylleraas Centre for Quantum Molecular Sciences), and from ERC-STG-2014 [grant number 639508].

Appendix. Proofs of main results

In this appendix all our proofs have been collected and are given in the same order as the corresponding results have appeared above. We will only present a partial proof of Theorem 2.1, as the proofs of points (1)  and (2)  are standard, and can be found in, e.g. Ref. [22].

Proof

Proof of Theorem 2.1

For point (3), we note that F is locally Lipschitz continuous as a consequence of $\mathcal{E}$ being smooth, which together with strong monotonicity makes points (1) and (2)  applicable. The remaining argument follows [23] closely (where the case $\mathcal{E}={\mathcal{E}}_{\mathrm{ECC}}$ was treated). First, by assumption of R, $F(Z_{\ast })=0$ and $F_d(Z_{\ast d})=0$ are equivalent to $\mathrm{\partial }\mathcal{E}(Z_{\ast })=0$ and $\mathrm{\partial }{\mathcal{E}}_d(Z_{\ast d})=0$, respectively (note that R commutes with $P_d$). Now, Taylor expanding $\mathcal{E}$ around $Z_{\ast }$ and evaluating at $Z_{\ast d}$ gives $\mathcal{E}(Z_{\ast d})-E_{\ast }=\frac{1}{2}⟨Z,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })Z⟩+\mathcal{O}(\Vert Z{\Vert }^3).$ By the smoothness of $\mathcal{E}$, there exists a constant $C^{\prime }$ such that $⟨Z,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })Z⟩\le C^{\prime }\Vert Z{\Vert }^2.$ Further, the fact that on U we can control the higher-order terms by the quadratic one, we have $|\mathcal{E}(Z_{\ast d})-E_{\ast }|\le C\Vert Z_{\ast d}-Z_{\ast }{\Vert }^2.$ Using Equation (8(8) $\Vert Z_{\ast d}-Z_{\ast }\Vert \le \frac{L}{\eta }\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d).$(8) ) gives the full statement in Equation (9(9) $\begin{array}{rl}|\mathcal{E}(Z_{\ast d})-E_{\ast }|& \le C\Vert Z_{\ast d}-Z_{\ast }{\Vert }^2\\ & \le C{\left(\frac{L}{\eta }\right)}^2\mathrm{dist}(Z_{\ast },{\mathcal{V}}_d\oplus {\mathcal{V}}_d{)}^2.\end{array}$(9) )

Proof

Proof of Theorem 3.1

With $F:\mathcal{V}\oplus \mathcal{V}\to {\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime }$ the flipped gradient and $X\in \mathcal{V}\oplus \mathcal{V}$, we obtain ${⟨X,F(Z)⟩}_{\mathcal{V}\oplus \mathcal{V},{\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime }}={⟨\mathrm{RX},\mathrm{\partial }\mathcal{E}(Z)⟩}_{\mathcal{V}\oplus \mathcal{V},{\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime }}.$ In the sequel, we omit the specification of the spaces in the dual pairing.

Let $Z_{\ast }\in \mathcal{V}\oplus \mathcal{V}$ be such that $\mathrm{\partial }\mathcal{E}(Z_{\ast })=0$, i.e. $F(Z_{\ast })=0$. Since F is smooth, local strong monotonicity of F is equivalent to $\mathrm{\partial }F(Z_{\ast })\in \mathcal{B}(\mathcal{V}\oplus \mathcal{V},{\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime })$ (the set of bounded linear operators) being coercive, i.e. there exists an ${\eta }_{\ast }>0$ such that $\mathrm{\Delta }(X)=⟨X,\mathrm{\partial F}(Z_{\ast })X⟩\ge {\eta }_{\ast }\Vert X{\Vert }_{\mathcal{V}\oplus \mathcal{V}}^2,$ where $\mathrm{\Delta }(X)$ is defined by the first equality. (The constant η in Equation (6(6) $⟨F(Z_1)-F(Z_2),Z_1-Z_2⟩\ge \eta \Vert Z_1-Z_2{\Vert }^2,$(6) ) approaches ${\eta }_{\ast }$ as the ball U given in Theorem 2.1 of the local strong monotonicity approaches a point.) To see this, we find an expression for $\mathrm{\Delta }(h)$ in terms of the energy map, (A1) $⟨X,F(Z_{\ast }+\mathrm{ϵX})⟩=⟨X,F(Z_{\ast })+ϵF^{\prime }(Z_{\ast };X)⟩+\mathcal{O}({ϵ}^2).$(A1) Here $F^{\prime }(Z_{\ast };X)$ is the directional derivative in the direction of X such that $\begin{array}{rl}⟨X,\mathrm{\partial F}(Z_{\ast })X⟩& =⟨X,F^{\prime }(Z_{\ast };X)⟩\\ & =\frac{\mathrm{d}}{\mathrm{d}ϵ}⟨\mathrm{RX},\mathrm{\partial }\mathcal{E}(Z_{\ast }+\mathrm{ϵX})⟩{|}_{ϵ=0}\\ & =⟨\mathrm{RX},{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })X⟩,\end{array}$ where ${\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })\in \mathcal{B}(\mathcal{V}\oplus \mathcal{V},{\mathcal{V}}^{\prime }\oplus {\mathcal{V}}^{\prime })$. By choosing ε small enough, Equation (A1(A1) $⟨X,F(Z_{\ast }+\mathrm{ϵX})⟩=⟨X,F(Z_{\ast })+ϵF^{\prime }(Z_{\ast };X)⟩+\mathcal{O}({ϵ}^2).$(A1) ) and strong monotonicity gives the coercivity claim. The logical implication also goes in the reverse direction. (This will be used below.)

Recall that ${\mathcal{E}}_{\theta }=\mathcal{E}\circ \theta$ and that $F_{\theta }$ denotes the flipped gradient of ${\mathcal{E}}_{\theta }$. We use that $F_{\theta }$ is locally strongly monotone at $W_{\ast }={\theta }^{-1}(Z_{\ast })$ if and only if ${\mathrm{\Delta }}_{\theta }$ is coercive, i.e. ${\mathrm{\Delta }}_{\theta }(X)=⟨X,\mathrm{\partial }{F}_{\theta }(W_{\ast })X⟩\ge {\eta }_{\theta }\Vert X{\Vert }^2,$ for some ${\eta }_{\theta }>0$ and all $X\in \mathcal{V}\oplus \mathcal{V}$. A straightforward application of the chain rule now gives $\begin{array}{rl}{\mathrm{\Delta }}_{\theta }(X)& =⟨M_{\ast }\mathrm{RX},{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩,\\ M_{\ast }& =\mathrm{\partial \theta }(W_{\ast })\in \mathcal{B}(\mathcal{V}\oplus \mathcal{V},\mathcal{V}\oplus \mathcal{V}).\end{array}$

We note that this is almost $\mathrm{\Delta }(M_{\ast }X)$. Indeed, (A2) $\begin{array}{rl}{\mathrm{\Delta }}_{\theta }(X)& =⟨R{M}_{\ast }X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩\\ & +⟨[M_{\ast },R]X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩\\ & =\mathrm{\Delta }(M_{\ast }X)+⟨(M_{\ast }R-R{M}_{\ast })X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩.\end{array}$(A2) In particular, if $M_{\ast }R=R_{\ast }M$ then the last term vanishes in the utmost right-hand side of Equation (A2(A2) $\begin{array}{rl}{\mathrm{\Delta }}_{\theta }(X)& =⟨R{M}_{\ast }X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩\\ & +⟨[M_{\ast },R]X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩\\ & =\mathrm{\Delta }(M_{\ast }X)+⟨(M_{\ast }R-R{M}_{\ast })X,{\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })(M_{\ast }X)⟩.\end{array}$(A2) ), and we obtain monotonicity of $F_{\theta }$ but with a modified constant.

In the case where $M_{\ast }R\ne R{M}_{\ast }$, we write $M_{\ast }=I+m_{\ast }$, and note that $M_{\ast }R-R{M}_{\ast }=m_{\ast }R-R{m}_{\ast }$. We obtain, (A3) $\begin{array}{rl}{\mathrm{\Delta }}_{\theta }(X)& \ge \eta \Vert M_{\ast }X{\Vert }^2-\Vert {\mathrm{\partial }}^2\mathcal{E}(Z_{\ast })\Vert \Vert M_{\ast }\Vert \Vert m_{\ast }\Vert \Vert X{\Vert }^2\\ & \ge [\eta \Vert M_{\ast }^{-1}{\Vert }^{-2}-C\Vert (1+\Vert m_{\ast }\Vert )\Vert m_{\ast }\Vert ]\Vert X{\Vert }^2.\end{array}$(A3) Here, we used that θ has a smooth inverse, implying $\Vert M_{\ast }X\Vert \ge \Vert M_{\ast }^{-1}{\Vert }^{-1}\Vert X\Vert$, and that $\Vert M_{\ast }\Vert \le I+\Vert m_{\ast }\Vert$.

Proof

Proof of Corollary 3.2

We consider the Jacobian of the coordinate map, which on block form reads (A4) $\mathrm{\partial }\theta (T^{\prime },{\mathrm{\Lambda }}^{\prime })=\left(\begin{array}{cc}\frac{\mathrm{\partial }T}{\mathrm{\partial }{T}^{\prime }}& \frac{\mathrm{\partial }T}{\mathrm{\partial }{\mathrm{\Lambda }}^{\prime }}\\ \frac{\mathrm{\partial \Lambda }}{\mathrm{\partial }{T}^{\prime }}& \frac{\mathrm{\partial \Lambda }}{\mathrm{\partial }{\mathrm{\Lambda }}^{\prime }}\end{array}\right).$(A4) For the map ${\theta }_n$ (see Equation (17(17) $(T,\mathrm{\Lambda })={\theta }_n(T^{\prime },{\mathrm{\Lambda }}^{\prime })=(T^{\prime },G_n({\mathrm{\Lambda }}^{\prime })),$(17) )), we first observe that by definition, $e^{\mathrm{\Lambda }}=e^{{\mathrm{\Lambda }}^{\prime }}+\mathcal{O}(\Vert {\mathrm{\Lambda }}^{\prime }{\Vert }^{n+1}),$ from which it follows, by taking the logarithm and expanding the logarithm around $e^{{\mathrm{\Lambda }}^{\prime }}$, which is a finite Taylor series, $\mathrm{\Lambda }={\mathrm{\Lambda }}^{\prime }+\mathcal{O}(\Vert {\mathrm{\Lambda }}^{\prime }{\Vert }^{n+1}).$ We obtain (A5) $\mathrm{\partial }{\theta }_n(T^{\prime },{\mathrm{\Lambda }}^{\prime })=\left(\begin{array}{cc}I& 0\\ 0& I+\mathcal{O}(\Vert {\mathrm{\Lambda }}^{\prime }{\Vert }^{n+1})\end{array}\right).$(A5) For the map ${\theta }_n\circ {\theta }_{\mathrm{C}-\mathrm{ECC}}$ we have, using the chain rule, (A6) $\mathrm{\partial }{\theta }_{\mathrm{C}-\mathrm{ECC}}(T^{\prime },{\mathrm{\Lambda }}^{\prime })=\left(\begin{array}{cc}A({\mathrm{\Lambda }}^{\prime })& {\mathrm{\partial }}_{{\mathrm{\Lambda }}^{\prime }}A({\mathrm{\Lambda }}^{\prime })T\\ 0& I+\mathcal{O}(\Vert {\mathrm{\Lambda }}^{\prime }{\Vert }^{n+1})\end{array}\right).$(A6) Here, $A({\mathrm{\Lambda }}^{\prime })$ is the linear transformation on $\mathcal{V}$ such that $S(T^{\prime };{\mathrm{\Lambda }}^{\prime })=A({\mathrm{\Lambda }}^{\prime })T^{\prime }$ (see Equation (12(12) $\mathrm{QT}{\varphi }_0=Q{e}^{-{{\mathrm{\Lambda }}^{\prime }}^{†}}{T}^{\prime }{\varphi }_0,Q=I-|{\varphi }_0⟩⟨{\varphi }_0|.$(12) )), i.e.$A({\mathrm{\Lambda }}^{\prime })$ can be expressed in terms of the matrix representation of $e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}$, (A7) $A({\mathrm{\Lambda }}^{\prime })T^{\prime }=\sum _{\mu ,\nu \in \mathcal{I}}{X}_{\mu }⟨{\varphi }_{\mu },e^{-({\mathrm{\Lambda }}^{\prime }{)}^{†}}{\varphi }_{\nu }⟩⟨{\varphi }_{\nu },T^{\prime }{\varphi }_0⟩.$(A7) We have $A({\mathrm{\Lambda }}^{\prime })=I+\mathcal{O}(\Vert {\mathrm{\Lambda }}^{\prime }\Vert )$, and ${\mathrm{\partial }}_{{\mathrm{\Lambda }}^{\prime }}A({\mathrm{\Lambda }}^{\prime })T^{\prime }=\mathcal{O}(\Vert T^{\prime }\Vert )$. For both maps, the Jacobian of the coordinate transformation at the critical point $W_{\ast }=(T_{\ast }^{\prime },{\mathrm{\Lambda }}_{\ast }^{\prime })$ becomes $M_{\ast }=I+m_{\ast }$ with $m_{\ast }=\mathcal{O}(\Vert W_{\ast }\Vert )$. Applying Theorem 3.1(2), the local strong monotonicity follows, and by Theorem 2.1, quasi-optimality of the truncated solutions and a quadratic error estimate.