CoupledCluster
Table of Contents
1. Brief description
This algorithm solves the CCSD amplitude equations using the Coulomb integrals and Hartree–Fock eigenenergies.
2. Algorithm call
A typical input file snippet to call the CoupledCluster algorithm is given below.
- name: CoupledCluster in: method: Ccsd integralsSliceSize: 100 slicedEigenEnergies: EigenEnergies coulombIntegrals: CoulombIntegrals slicedCoulombVertex: CoulombVertex maxIterations: 20 energyConvergence: 1.0E-4 amplitudesConvergence: 1.0E-4 mixer: type: DiisMixer maxResidua: 5 #or for instance a linear mixer #mixer: # type: LinearMixer # ratio: 1.0 out: amplitudes: Amplitudes
3. Algorithm input
| Keyword | Value |
|---|---|
method |
Type of coupled cluster algorithm |
linearized |
Use linearized or fully non-linear coupled cluster |
integralsSliceSize |
Integral slice size |
slicedEigenEnergies |
Sliced one-electron energies |
coulombIntegrals |
Coulomb Integrals |
slicedCoulombVertex |
Sliced Coulomb vertex |
maxIterations |
Maximum number of iterations |
energyConvergence |
Energy convergence threshold |
amplitudesConvergence |
Amplitude convergence threshold |
mixer |
Mixer type and parameters |
initialAmplitudes |
initial Amplitudes to begin the iterations |
3.1. method
method specifies which approximation of the coupled-cluster
Amplitudes should be used. See the
Amplitudes for the specification of
each approximation.
The currently available approximations are:
Ccsd- Coupled Cluster singles and doubles. See (Bartlett and Musiał 2007) and references therin.
CcsdReference- Canonical and unoptimized version of the
Ccsdalgorithm. This algorithm has the same inputs asCcsdwith the exception that the arguments integralSliceSize and slicedCoulombVertex and are ignored. Drccd- Unoptimized implementation of Direct-Ring Coupled Cluster singles and doubles giving the Random Phase Approximation (RPA) with the Second Order Screened Exchange (SOSEX) correction (Freeman 1977; Grüneis et al. 2009).
3.2. linearzied
linearized A non-zero value indicates that the linearized
coupled-cluster amplitude equations should be solved, rather
than the fully non-linear ones.
The default is linearized: 0.
Note: currently only method: Drccd supports linearized
amplitude equations.
3.3. integralsSliceSize
integralsSliceSize controls the slice size of the \(V_{cd}^{ab}\) integrals, which are computed on-the-fly to
reduce the memory footprint. The integer value specified for integralsSliceSize refers to the dimension size
used for the \(a\) and \(b\) index. We recommend setting integralsSliceSize: 100 to balance
computational efficency with memory usage. Smaller/larger values reduce/increase the memory footprint.
3.4. maxIterations
maxIterations controls the maximum number of iterations allowed to solve the \(t_{ij}^{ab}\) and \(t_i^a\) amplitude equations.
If convergence of the energy and residual vectors within the specified thresholds is achieved using fewer iterations
than maxIterations, the algorithm was successfull and will stop.
If maxIterations is reached without achieving energyConvergence and amplitudesConvergence , the algorithm was not
successful and will stop. We recommend to set maxIterations: 20, which is ususally enough to achieve reasonable convergence
thresholds.
3.5. energyConvergence
energyConvergence specifies the convergence threshold for the correlation energy.
If energyConvergence and amplitudesConvergence is achieved, the iterative solution was successful and the algorithm will stop.
3.6. amplitudesConvergence
amplitudesConvergence specifies the convergence threshold for the residual vector of the singles and doubles amplitude equations.
If energyConvergence and amplitudesConvergence is achieved, the iterative solution was successful and the algorithm will stop.
3.7. initialAmplitudes
initialAmplitudes specifies amplitudes from a precedent calculations. This allows restarting a self consistent calculation from a former checkpoint.
3.8. mixer
mixer specifies mixer-specific parameters used to solve the amplitude equations iteratively.
The direct inversion iterative subspace (Diis) algorithm is the default algorithm used to mix amplitude guesses from previous
iterations to improve the guess for the next iteration.
mixer: type can currently be set to *DiisMixer or *LinearMixer.
3.8.1. DiisMixer
For mixer: type: DissMixer, it is possible to specify mixer: maxResidua, which controls the number of residual
vectors used by the Diis mixer.
We recommend to set the maximum number of residual vectors used in the DiisMixer to 5.
More residual vectors result in a larger memory footprint.
- Example
mixer: type: DiisMixer maxResidua: 5
3.8.2. LinearMixer
For mixer: type: LinearMixer it is possible to specify mixer: ratio, which controls the mixing ratio used
to update the new guess of the amplitudes based on estimates from the previous and current iteration.
The LinearMixer exhibits the smallest possible memory footprint from all currently available mixers,
keeping two sets of amplitudes in memory at once only.
We recommend to set mixer: ratio: 1.0. Smaller ratios result in a slower but perhaps more stable convergence.
- Example
mixer: type: LinearMixer ratio: 1.0
4. Algorithm output
| Keyword | Value |
|---|---|
amplitudes |
Amplitudes |
The output of the CoupledCluster algorithm includes energy and amplitudes . The amplitudes output contains
the converged singles and doubles amplitude tensors. The amplitudes can be used as input for algorithms
that estimate the finite simulation cell size error (FiniteSizeCorrection)
and the basis set incompleteness errors (BasisSetCorrection).
4.1. Sample stdout
Below an example standard output stream is shown for a successful CoupledCluster algorithm run.
step: 6, CoupledCluster Using method Ccsd. integralsSliceSize: 100 Using mixer DiisMixer. maxResidua: 5 Maximum number of iterations: 30 Unless reaching energy convergence dE: 1e-05 and amplitudes convergence dR: 1e-05 Iter Energy dE dR time GF/s/core 1 -2.43605043e+01 -2.4361e+01 4.3924e-01 0.5 1.0 2 -2.47577534e+01 -3.9725e-01 7.4733e-02 0.8 4.8 3 -2.53776918e+01 -6.1994e-01 1.8674e-02 0.7 5.4 4 -2.54455925e+01 -6.7901e-02 6.4132e-03 0.7 5.5 5 -2.54445080e+01 1.0845e-03 2.2120e-03 0.7 5.4 6 -2.54458312e+01 -1.3232e-03 1.0304e-03 0.7 5.4 7 -2.54448941e+01 9.3705e-04 5.0727e-04 0.7 5.4 8 -2.54452894e+01 -3.9521e-04 1.9694e-04 0.7 5.4 9 -2.54454262e+01 -1.3682e-04 7.7180e-05 0.7 5.4 10 -2.54455328e+01 -1.0663e-04 3.0247e-05 0.7 5.4 11 -2.54455929e+01 -6.0110e-05 1.1758e-05 0.7 5.4 12 -2.54456151e+01 -2.2218e-05 5.0053e-06 0.7 5.4 13 -2.54456249e+01 -9.7454e-06 2.2689e-06 0.7 5.4 Ccsd correlation energy: -25.4456248862 2nd-order correlation energy: -24.3605043096 realtime 9.189542891 s --
5. Sample yaml output
Below an example yaml output stream is shown for a successful CoupledCluster algorithm run.
name: CoupledCluster out: amplitudes: 0x26e4758 convergenceReached: 1 energy: correlation: -25.445624886202758 direct: -38.822491455744313 exchange: 13.376866569541555 secondOrder: -24.360504309639897 unit: 0.036749322175638782 realtime: 9.189542891
6. Computational complexity
This section explains computational or memory footprints for the various methods implemented in CoupledCluster (see method).
CcsdThe computational bottle neck of a typical Ccsd calculation originates from the following contraction, which is part of the doubles amplitude equations: \(V_{cd}^{ab} t_{ij}^{cd}\). The computational cost for evaluating this expression scales as \(\mathcal{O}(N_{\rm o}^2 N_{\rm v}^4)\). To avoid a memory footprint of \(\mathcal{O}(N_{\rm v}^4)\) in storing \(V_{cd}^{ab}\), slices \(V_{cd}^{xy}\) are computed on-the-fly and used in the contraction, reducing the corresponding memory footprint to \(\mathcal{O}(N_{\rm v}^2 N_{\rm s}^2)\), where \(N_{\rm s}\) is controlled using the
integralsSliceSizekeyword.We note that required storage of a set of doubles amplitudes adds substantially to the memory footprint in Ccsd calculations. The Diis algorithm requires the storage of both doubles residua and amplitudes
maxResiduatimes. We recommend to choose the type of mixer and its parameters carefully to reduce the memory footprint if necessary.Drccd- The computational complexity is \(\mathcal O(N_\mathrm{o}^3N_\mathrm{v}^3)\). The implementation is not optimized for large systems. The memory requirement scales as \(\mathcal O(N_\mathrm{o}^2 N_\mathrm{v}^2)\)
7. Theory
Coupled-cluster employs the exponential Ansatz for the correlated wave function
\begin{equation} | \Psi \rangle = e^{\hat T} | \Phi \rangle, \end{equation}where \(|\Phi\rangle\) denotes the single Hartree–Fock slater determinant. The cluster operator \(\hat T = \hat T_1 + \hat T_2 + \ldots\) is expanded in increasing number of excitation levels. The single and double exciation parts of the cluster operator are given by
\begin{eqnarray} \hat T_1 = \sum_{ai} t^a_i \hat\tau^a_i, \\ \hat T_2 = \sum_{abij} t^{ab}_{ij} \hat\tau^{ab}_{ij}, \end{eqnarray}
where \(\hat \tau^{a\ldots}_{i\ldots} = \hat c^\dagger_a\ldots \hat c_i\ldots\)
denotes the exciation operator.
The coefficients \(t^a_i\), \(t^{ab}_{ij}\), \(\ldots\) are called coupled-cluster
Amplitudes.
They are found by projecting the stationary Schrödinger equation for the
coupled-cluster Ansatz \(E|\Psi\rangle = \hat H|\Psi\rangle\)
on \(\langle \Phi|e^{-\hat T}\).
The resulting over-determined set of equations is truncated at the
same level as the cluster-operator \(\hat T\), giving as many equations
as there are coefficients in \(\hat T\).
The Amplitudes are generated by the CoupledCluster
algorithm by solving the amplitude equation of the employed
coupled-cluster method, described below:
7.1. Coupled-Cluster Singles Doubles (Ccsd)
The cluster-operator is truncated containing only single and double excitations. The projections on the singly and doubly excited slater-determinants \(\langle \Phi^a_i|\) and \(\langle \Phi^{ab}_{ij}|\) give
\begin{align} \big\langle \Phi^{a}_{i} \big| e^{-\hat T} \hat H e^{\hat T} \big| \Phi \big\rangle &= 0 \quad \forall ai, \\ \label{eqn:t2} \big\langle \Phi^{ab}_{ij} \big| e^{-\hat T} \hat H e^{\hat T} \big| \Phi \big\rangle &= 0 \quad \forall abij. \end{align}The above equations expand to a finite set of a few dozens of, in general, non-linear algebraic contractions of \(t^a_i\) and \(t^{ab}_{ij}\) with \(\varepsilon_p\) and \(V^{pq}_{sr}\). See (Shavitt and Bartlett 2009) for a text-book introduction of the original works of (Coester and Kümmel 1960) and (Čížek 1969).
7.2. Direct-Ring Couple-Cluster Doubles (Drccd)
This method uses only the double excitation part of the cluster-operator. From the full doubles amplitude equations in Eq. (\ref{eqn:t2}) only those contractions are considered where pairs of particle and hole contractions originate from a common vertex and terminate in another common vertex, forming direct-rings in the diagrammatic notation. Only five terms remain in a canonical calculation and they read
\begin{equation} \Delta^{ab}_{ij} t^{ab}_{ij} + V^{ab}_{ij} + V^{kb}_{cj} t^{ac}_{ik} + V^{al}_{id} t^{db}_{lj} + V^{kl}_{cd} t^{ac}_{ik} t^{db}_{lj} = 0 \quad \forall abij, \end{equation}
with
\(\Delta^{ab}_{ij} = \varepsilon_a+\varepsilon_b-\varepsilon_i-\varepsilon_j\)
and where a sum over all indices that occurr only on the left-hand-side
is implied. These terms are the dominant terms of coupled-cluster singles
doubles in the hight-density limit (Gell-Mann and Brueckner 1957).
The direct correlation contribution of the Drccd method
\(\sum_{abij} 2t^{ab}_{ij}V^{ij}_{ab}\) gives the Random Phase Approximation (RPA),
the remaining exchange contribution is often termed Second Order
Screened Exchange (SOSEX) correction
(Freeman 1977; Grüneis et al. 2009).
See (Furche 2008; Chen, Agee, and Furche 2018) for a review on the RPA and its corrections,
as well as (Macke 1950; Pines and Bohm 1952) for the
original work on the RPA.
We recommend the following review article and references therein to get started with coupled-cluster theory (Bartlett and Musiał 2007) .