14. cc — Coupled cluster

The cc module implements the coupled cluster (CC) model to compute energies, analytical nuclear gradients, density matrices, excited states, and relevant properties.

To compute the CC energy, one first needs to perform a mean-field calculation using the mean-field module scf. The mean-field object defines the Hamiltonian and the problem size, which are used to initialize the CC object:

from pyscf import gto, scf, cc
mol = gto.M(atom='H 0 0 0; F 0 0 1', basis='ccpvdz')
mf = scf.RHF(mol).run()
mycc = cc.CCSD(mf)
mycc.kernel()

Unrelaxed density matrices are evaluated in the MO basis:

dm1 = mycc.make_rdm1()
dm2 = mycc.make_rdm2()

The CCSD(T) energy can be obtained by:

from pyscf.cc import ccsd_t
print(ccsd_t.kernel(mycc, mycc.ao2mo())[0])

Gradients are available:

from pyscf.cc import ccsd_grad
from pyscf import grad
grad_e = ccsd_grad.kernel(mycc)
grad_n = grad.grad_nuc(mol)
grad = grad_e + grad_nuc

Excited states can be calculated with ionization potential (IP), electron affinity (EA), and electronic excitation (EE) equation-of-motion (EOM) CCSD:

mycc = cc.RCCSD(mf)
mycc.kernel()
e_ip, c_ip = mycc.ipccsd(nroots=1)
e_ea, c_ea = mycc.eaccsd(nroots=1)
e_ee, c_ee = mycc.eeccsd(nroots=1)

mycc = cc.UCCSD(mf)
mycc.kernel()
e_ip, c_ip = mycc.ipccsd(nroots=1)
e_ea, c_ea = mycc.eaccsd(nroots=1)
e_ee, c_ee = mycc.eeccsd(nroots=1)

All CC methods have two implementations. One is simple and highly readable (suffixed by _slow in the filename) and the other is extensively optimized for computational efficiency. All code in the _slow versions is structured as close as possible to the formulas documented in the literature. Pure Python/numpy data structures and functions are used so that explicit memory management is avoided. It is easy to make modifications or develop new methods based on the slow implementations.

The computationally efficient (outcore) version is the default implementation for the CC module. In this implementation, the CPU usage, memory footprint, memory efficiency, and IO overhead are carefully considered. To keep a small memory footprint, most integral tensors are stored on disk. IO is one of the main bottlenecks in this implementation. Two techniques are used to reduce the IO overhead. One is the asynchronized IO to overlap the computation and reading/writing of the 4-index tensors. The other is AO-driven for the contraction of T2 and (vv|vv) integrals in CCSD and CCSD-lambda functions. These techniques allow the CC module to efficiently handle medium-sized systems. In a test system with 25 occupied orbitals and 1500 virtual orbitals, each CCSD iteration takes about 2.5 hours. The program does not automatically switch to AO-driven CCSD for large systems. The user must manually set the direct attribute to enable an AO-driven CCSD calculation:

mycc = cc.CCSD(mf)
mycc.direct = True
mycc.kernel()

Some of the CC methods have an efficient incore implementation, where all tensors are held in memory. The incore implementation reduces the IO overhead and optimizes certain formulas to gain the best FLOPS. It is about 30% faster than the outcore implementation. Depending on the available memory, the incore code can be used for systems with up to approximately 250 orbitals.

Point group symmetry is not considered in the CCSD programs, but it is used in the CCSD(T) code to gain the best performance.

Arbitrary frozen orbitals (not limited to frozen core) are supported by the CCSD, CCSD(T), density matrices, and EOM-CCSD modules, but not in the analytical CCSD gradient module.

14.1. Examples

This section documents some examples about how to effectively use the CCSD module, and how to incorporate the CCSD solver with other PySCF functions to perform advanced simulations. For a complete list of CC examples, see pyscf/examples/cc.

14.1.1. A general solver for customized Hamiltonian

The CC module is not limited to molecular systems. The program is implemented as a general solver for arbitrary Hamiltonians. It allows users to overwrite the default molecular Hamiltonian with their own effective Hamiltonians. In this example, we create a Hubbard model and feed its Hamiltonian to the CCSD module.

#!/usr/bin/env python

'''
Six-site 1D U/t=2 Hubbard-like model system with PBC at half filling.
The model is gapped at the mean-field level
'''

import numpy
from pyscf import gto, scf, ao2mo, cc

mol = gto.M(verbose=4)
n = 6
mol.nelectron = n
# Setting incore_anyway=True to ensure the customized Hamiltonian (the _eri
# attribute) to be used in the post-HF calculations.  Without this parameter,
# some post-HF method (particularly in the MO integral transformation) may
# ignore the customized Hamiltonian if memory is not enough.
mol.incore_anyway = True

h1 = numpy.zeros((n,n))
for i in range(n-1):
    h1[i,i+1] = h1[i+1,i] = -1.0
h1[n-1,0] = h1[0,n-1] = -1.0
eri = numpy.zeros((n,n,n,n))
for i in range(n):
    eri[i,i,i,i] = 2.0

mf = scf.RHF(mol)
mf.get_hcore = lambda *args: h1
mf.get_ovlp = lambda *args: numpy.eye(n)
mf._eri = ao2mo.restore(8, eri, n)
mf.kernel()


# In PySCF, the customized Hamiltonian needs to be created once in mf object.
# The Hamiltonian will be used everywhere whenever possible.  Here, the model
# Hamiltonian is passed to CCSD object via the mf object.

mycc = cc.RCCSD(mf)
mycc.kernel()
e,v = mycc.ipccsd(nroots=3)
print(e)

14.1.2. Using CCSD as CASCI active space solver

CCSD program can be wrapped as a Full CI solver, which can be combined with the CASCI solver to approximate the multi-configuration calculation.

#!/usr/bin/env python

'''
Using the CCSD method as the active space solver to compute an approximate
CASCI energy.

A wrapper is required to adapt the CCSD solver to CASCI fcisolver interface.
Inside the wrapper function, the CCSD code is the same as the example
40-ccsd_with_given_hamiltonian.py
'''

import numpy
from pyscf import gto, scf, cc, ao2mo, mcscf

class AsFCISolver(object):
    def __init__(self):
        self.mycc = None

    def kernel(self, h1, h2, norb, nelec, ci0=None, ecore=0, **kwargs):
        fakemol = gto.M(verbose=0)
        nelec = numpy.sum(nelec)
        fakemol.nelectron = nelec
        fake_hf = scf.RHF(fakemol)
        fake_hf._eri = ao2mo.restore(8, h2, norb)
        fake_hf.get_hcore = lambda *args: h1
        fake_hf.get_ovlp = lambda *args: numpy.eye(norb)
        fake_hf.kernel()
        self.mycc = cc.CCSD(fake_hf)
        self.eris = self.mycc.ao2mo()
        e_corr, t1, t2 = self.mycc.kernel(eris=self.eris)
        l1, l2 = self.mycc.solve_lambda(t1, t2, eris=self.eris)
        e = fake_hf.e_tot + e_corr
        return e+ecore, [t1,t2,l1,l2]

    def make_rdm1(self, fake_ci, norb, nelec):
        mo = self.mycc.mo_coeff
        t1, t2, l1, l2 = fake_ci
        dm1 = reduce(numpy.dot, (mo, self.mycc.make_rdm1(t1, t2, l1, l2), mo.T))
        return dm1

    def make_rdm12(self, fake_ci, norb, nelec):
        mo = self.mycc.mo_coeff
        nmo = mo.shape[1]
        t1, t2, l1, l2 = fake_ci
        dm2 = self.mycc.make_rdm2(t1, t2, l1, l2)
        dm2 = numpy.dot(mo, dm2.reshape(nmo,-1))
        dm2 = numpy.dot(dm2.reshape(-1,nmo), mo.T)
        dm2 = dm2.reshape([nmo]*4).transpose(2,3,0,1)
        dm2 = numpy.dot(mo, dm2.reshape(nmo,-1))
        dm2 = numpy.dot(dm2.reshape(-1,nmo), mo.T)
        dm2 = dm2.reshape([nmo]*4)
        return self.make_rdm1(fake_ci, norb, nelec), dm2

    def spin_square(self, fake_ci, norb, nelec):
        return 0, 1

mol = gto.M(atom = 'H 0 0 0; F 0 0 1.2',
            basis = 'ccpvdz',
            verbose = 4)
mf = scf.RHF(mol).run()
norb = mf.mo_coeff.shape[1]
nelec = mol.nelectron
mc = mcscf.CASCI(mf, norb, nelec)
mc.fcisolver = AsFCISolver()
mc.kernel()

14.1.3. Gamma point CCSD with Periodic boundary condition

Integrals in Gamma point of periodic Hartree-Fock calculation are all real. You can feed the integrals into any pyscf molecular module using the same operations as the above example. However, the interface between PBC code and molecular code are more compatible. You can treat the crystal object and the molecule object in the same manner. In this example, you can pass the PBC mean field method to CC module to have the gamma point CCSD correlation.

#!/usr/bin/env python

'''
Gamma point post-HF calculation needs only real integrals.
Methods implemented in finite-size system can be directly used here without
any modification.
'''

import numpy
from pyscf.pbc import gto, scf

cell = gto.M(
    a = numpy.eye(3)*3.5668,
    atom = '''C     0.      0.      0.    
              C     0.8917  0.8917  0.8917
              C     1.7834  1.7834  0.    
              C     2.6751  2.6751  0.8917
              C     1.7834  0.      1.7834
              C     2.6751  0.8917  2.6751
              C     0.      1.7834  1.7834
              C     0.8917  2.6751  2.6751''',
    basis = '6-31g',
    verbose = 4,
)

mf = scf.RHF(cell).density_fit()
mf.with_df.mesh = [10]*3
mf.kernel()

#
# Import CC, TDDFT moduel from the molecular implementations
#
from pyscf import cc, tddft
mycc = cc.CCSD(mf)
mycc.kernel()

mytd = tddft.TDHF(mf)
mytd.nstates = 5
mytd.kernel()

14.1.4. CCSD with truncated MOs to avoid linear dependency

It is common to have linear dependence when one wants to systematically enlarge the AO basis set to approach complete basis set limit. The numerical instability usually has noticeable effects on the CCSD convergence. An effective way to remove this negative effects is to truncate the AO sets and allow the MO orbitals being less than AO functions.

#!/usr/bin/env python

'''
:func:`scf.addons.remove_linear_dep_` discards the small eigenvalues of overlap
matrix.  This reduces the number of MOs from 50 to 49.  The problem size of
the following CCSD method is 49.
'''

from pyscf import gto, scf, cc
mol = gto.Mole()
mol.atom = [('H', 0, 0, .5*i) for i in range(20)]
mol.basis = 'ccpvdz'
mol.verbose = 4
mol.build()
mf = scf.RHF(mol).run()
mycc = cc.CCSD(mf).run()

mf = scf.addons.remove_linear_dep_(mf).run()
mycc = cc.CCSD(mf).run()

14.1.5. Response and un-relaxed CCSD density matrix

CCSD has two kinds of one-particle density matrices. The (second order) un-relaxed density matrix and the (relaxed) response density matrix. The CCSD.make_rdm1() function computes the un-relaxed density matrix which is associated to the regular CCSD energy formula. The response density is mainly used to compute the first order response quantities eg the analytical nuclear gradients. It is not recommended to use the response density matrix for population analysis.

#!/usr/bin/env python
#
# Author: Qiming Sun <osirpt.sun@gmail.com>
#

'''
CCSD density matrix
'''

import numpy
from pyscf import gto, scf, cc, ao2mo

mol = gto.M(
    atom = 'H 0 0 0; F 0 0 1.1',
    basis = 'ccpvdz')
mf = scf.RHF(mol).run()
mycc = cc.CCSD(mf).run()

#
# CCSD density matrix in MO basis
#
dm1 = mycc.make_rdm1()
dm2 = mycc.make_rdm2()

#
# CCSD energy based on density matrices
#
h1 = numpy.einsum('pi,pq,qj->ij', mf.mo_coeff.conj(), mf.get_hcore(), mf.mo_coeff)
nmo = mf.mo_coeff.shape[1]
eri = ao2mo.kernel(mol, mf.mo_coeff, compact=False).reshape([nmo]*4)
E = numpy.einsum('pq,qp', h1, dm1)
# Note dm2 is transposed to simplify its contraction to integrals
E+= numpy.einsum('pqrs,pqrs', eri, dm2) * .5
E+= mol.energy_nuc()
print('E(CCSD) = %s, reference %s' % (E, mycc.e_tot))


# When plotting CCSD density on the grids, CCSD density matrices need to be
# first transformed to AO basis.
dm1_ao = numpy.einsum('pi,ij,qj->pq', mf.mo_coeff, dm1, mf.mo_coeff.conj())

from pyscf.tools import cubegen
cubegen.density(mol, 'rho_ccsd.cube', dm1_ao)

14.1.6. Reusing integrals in CCSD and relevant calculations

By default the CCSD solver and the relevant CCSD lambda solver, CCSD(T), CCSD gradients program generate MO integrals in their own runtime. But in most scenario, the same MO integrals can be generated once and reused in the four modules. To remove the overhead of recomputing MO integrals, the three module support user to feed MO integrals.

#!/usr/bin/env python
#
# Author: Qiming Sun <osirpt.sun@gmail.com>
#

'''
To avoid recomputing AO to MO integral transformation, integrals for CCSD,
CCSD(T), CCSD lambda equation etc can be reused.
'''

from pyscf import gto, scf, cc

mol = gto.M(verbose = 4,
    atom = 'H 0 0 0; F 0 0 1.1',
    basis = 'ccpvdz')

mf = scf.RHF(mol).run()
mycc = cc.CCSD(mf)

#
# CCSD module allows you feed MO integrals
#
eris = mycc.ao2mo()
mycc.kernel(eris=eris)

#
# The same MO integrals can be used in CCSD lambda equation
#
mycc.solve_lambda(eris=eris)

#
# CCSD(T) module requires the same integrals used by CCSD module
#
from pyscf.cc import ccsd_t
ccsd_t.kernel(mycc, eris=eris)

#
# CCSD gradients need regular MO integrals to solve the relaxed 1-particle
# density matrix
#
from pyscf.cc import ccsd_grad
grad_e = ccsd_grad.kernel(mycc, eris=eris)  # The electronic part only

14.1.7. Interfering CCSD-DIIS

14.1.8. Restart CCSD

14.2. Program reference

14.2.1. cc.ccsd module and CCSD class

The pyscf.cc.ccsd.CCSD class is the object to hold the restricted CCSD environment attributes and results. The environment attributes are the parameters to control the runtime behavior of the CCSD module, e.g. the convergence criteria, DIIS parameters, and so on. After the ground state CCSD calculation, correlation energy, T1 and T2 amplitudes are stored in the CCSD object. This class supports the calculation of CCSD 1- and 2-particle density matrices.

class pyscf.cc.ccsd.CCSD(mf, frozen=0, mo_coeff=None, mo_occ=None)[source]

restricted CCSD

Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory
conv_tol
: float
converge threshold. Default is 1e-7.
conv_tol_normt
: float
converge threshold for norm(t1,t2). Default is 1e-5.
max_cycle
: int
max number of iterations. Default is 50.
diis_space
: int
DIIS space size. Default is 6.
diis_start_cycle
: int
The step to start DIIS. Default is 0.
iterative_damping
: float
The self consistent damping parameter.
direct
: bool
AO-direct CCSD. Default is False.
async_io
: bool
Allow for asynchronous function execution. Default is True.
incore_complete
: bool
Avoid all I/O (also for DIIS). Default is False.
level_shift
: float
A shift on virtual orbital energies to stablize the CCSD iteration
frozen
: int or list

If integer is given, the inner-most orbitals are frozen from CC amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in CC calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> mycc = cc.CCSD(mf).set(frozen = 2).run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> mycc.set(frozen = [0,1,16,17,18]).run()

Saved results

converged
: bool
CCSD converged or not
e_corr
: float
CCSD correlation correction
e_tot
: float
Total CCSD energy (HF + correlation)
t1, t2 :
T amplitudes t1[i,a], t2[i,j,a,b] (i,j in occ, a,b in virt)
l1, l2 :
Lambda amplitudes l1[i,a], l2[i,j,a,b] (i,j in occ, a,b in virt)

RCCSD for real integrals 8-fold permutation symmetry has been used (ij|kl) = (ji|kl) = (kl|ij) = ...

pyscf.cc.ccsd.CC

alias of CCSD

class pyscf.cc.ccsd.CCSD(mf, frozen=0, mo_coeff=None, mo_occ=None)[source]

restricted CCSD

Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory
conv_tol
: float
converge threshold. Default is 1e-7.
conv_tol_normt
: float
converge threshold for norm(t1,t2). Default is 1e-5.
max_cycle
: int
max number of iterations. Default is 50.
diis_space
: int
DIIS space size. Default is 6.
diis_start_cycle
: int
The step to start DIIS. Default is 0.
iterative_damping
: float
The self consistent damping parameter.
direct
: bool
AO-direct CCSD. Default is False.
async_io
: bool
Allow for asynchronous function execution. Default is True.
incore_complete
: bool
Avoid all I/O (also for DIIS). Default is False.
level_shift
: float
A shift on virtual orbital energies to stablize the CCSD iteration
frozen
: int or list

If integer is given, the inner-most orbitals are frozen from CC amplitudes. Given the orbital indices (0-based) in a list, both occupied and virtual orbitals can be frozen in CC calculation.

>>> mol = gto.M(atom = 'H 0 0 0; F 0 0 1.1', basis = 'ccpvdz')
>>> mf = scf.RHF(mol).run()
>>> # freeze 2 core orbitals
>>> mycc = cc.CCSD(mf).set(frozen = 2).run()
>>> # freeze 2 core orbitals and 3 high lying unoccupied orbitals
>>> mycc.set(frozen = [0,1,16,17,18]).run()

Saved results

converged
: bool
CCSD converged or not
e_corr
: float
CCSD correlation correction
e_tot
: float
Total CCSD energy (HF + correlation)
t1, t2 :
T amplitudes t1[i,a], t2[i,j,a,b] (i,j in occ, a,b in virt)
l1, l2 :
Lambda amplitudes l1[i,a], l2[i,j,a,b] (i,j in occ, a,b in virt)
as_scanner(cc)

Generating a scanner/solver for CCSD PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CCSD energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CCSD and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, cc
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> cc_scanner = cc.CCSD(scf.RHF(mol)).as_scanner()
>>> e_tot = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
energy(mycc, t1=None, t2=None, eris=None)

CCSD correlation energy

get_frozen_mask(mp)

Get boolean mask for the restricted reference orbitals.

In the returned boolean (mask) array of frozen orbital indices, the element is False if it corresonds to the frozen orbital.

make_rdm1(t1=None, t2=None, l1=None, l2=None, ao_repr=False)[source]

Un-relaxed 1-particle density matrix in MO space

make_rdm2(t1=None, t2=None, l1=None, l2=None)[source]

2-particle density matrix in MO space. The density matrix is stored as

dm2[p,r,q,s] = <p^+ q^+ s r>

restore_from_diis_(mycc, diis_file, inplace=True)

Reuse an existed DIIS object in the CCSD calculation.

The CCSD amplitudes will be restored from the DIIS object to generate t1 and t2 amplitudes. The t1/t2 amplitudes of the CCSD object will be overwritten by the generated t1 and t2 amplitudes. The amplitudes vector and error vector will be reused in the CCSD calculation.

pyscf.cc.ccsd.RCCSD

alias of CCSD

pyscf.cc.ccsd.as_scanner(cc)[source]

Generating a scanner/solver for CCSD PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CCSD energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CCSD and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, cc
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> cc_scanner = cc.CCSD(scf.RHF(mol)).as_scanner()
>>> e_tot = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))
pyscf.cc.ccsd.energy(mycc, t1=None, t2=None, eris=None)[source]

CCSD correlation energy

pyscf.cc.ccsd.restore_from_diis_(mycc, diis_file, inplace=True)[source]

Reuse an existed DIIS object in the CCSD calculation.

The CCSD amplitudes will be restored from the DIIS object to generate t1 and t2 amplitudes. The t1/t2 amplitudes of the CCSD object will be overwritten by the generated t1 and t2 amplitudes. The amplitudes vector and error vector will be reused in the CCSD calculation.

14.2.2. cc.rccsd and RCCSD class

pyscf.cc.rccsd.RCCSD is also a class for restricted CCSD calculations, but different to the pyscf.cc.ccsd.CCSD class. It uses different formula to compute the ground state CCSD solution. Although slower than the implmentation in the pyscf.cc.ccsd.CCSD class, it supports the system with complex integrals. Another difference is that this class supports EOM-CCSD methods, including EOM-IP-CCSD, EOM-EA-CCSD, EOM-EE-CCSD, EOM-SF-CCSD.

class pyscf.cc.rccsd.RCCSD(mf, frozen=0, mo_coeff=None, mo_occ=None)[source]

restricted CCSD with IP-EOM, EA-EOM, EE-EOM, and SF-EOM capabilities

Ground-state CCSD is performed in optimized ccsd.CCSD and EOM is performed here.

Restricted CCSD implementation which supports both real and complex integrals. The 4-index integrals are saved on disk entirely (without using any symmetry). This code is slower than the pyscf.cc.ccsd implementation.

Note MO integrals are treated in chemist’s notation

Ref: Hirata et al., J. Chem. Phys. 120, 2581 (2004)

class pyscf.cc.rccsd.RCCSD(mf, frozen=0, mo_coeff=None, mo_occ=None)[source]

restricted CCSD with IP-EOM, EA-EOM, EE-EOM, and SF-EOM capabilities

Ground-state CCSD is performed in optimized ccsd.CCSD and EOM is performed here.

ccsd(t1=None, t2=None, eris=None, mbpt2=False)[source]

Ground-state CCSD.

Kwargs:
mbpt2
: bool
Use one-shot MBPT2 approximation to CCSD.
energy(cc, t1=None, t2=None, eris=None)

RCCSD correlation energy

pyscf.cc.rccsd.energy(cc, t1=None, t2=None, eris=None)[source]

RCCSD correlation energy

14.2.3. cc.uccsd and UCCSD class

pyscf.cc.uccsd.UCCSD class supports the CCSD calculation based on UHF wavefunction as well as the ROHF wavefunction. Besides the ground state UCCSD calculation, UCCSD lambda equation, 1-particle and 2-particle density matrices, EOM-IP-CCSD, EOM-EA-CCSD, EOM-EE-CCSD are all available in this class. Note this class does not support complex integrals.

class pyscf.cc.uccsd.UCCSD(mf, frozen=0, mo_coeff=None, mo_occ=None)[source]

UCCSD with spatial integrals

pyscf.cc.uccsd.energy(cc, t1=None, t2=None, eris=None)[source]

UCCSD correlation energy

14.2.4. cc.addons

Helper functions for CCSD, RCCSD and UCCSD modules are implemented in cc.addons

pyscf.cc.addons.spatial2spin(tx, orbspin=None)[source]

Convert T1/T2 of spatial orbital representation to T1/T2 of spin-orbital representation

pyscf.cc.addons.spatial2spinorb(tx, orbspin=None)

Convert T1/T2 of spatial orbital representation to T1/T2 of spin-orbital representation

14.2.5. CCSD(T)

RHF-CCSD(T) for real integrals

14.2.6. CCSD gradients

pyscf.cc.ccsd_grad.as_scanner(cc)

Generating a scanner/solver for CCSD PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CCSD energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters assigned in the CCSD and the underlying SCF objects (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, cc
>>> mol = gto.M(atom='H 0 0 0; F 0 0 1')
>>> cc_scanner = cc.CCSD(scf.RHF(mol)).as_scanner()
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.1'))
>>> e_tot, grad = cc_scanner(gto.M(atom='H 0 0 0; F 0 0 1.5'))