Skip to content

Latest commit

 

History

History
406 lines (309 loc) · 14.4 KB

df.rst

File metadata and controls

406 lines (309 loc) · 14.4 KB

:mod:`df` --- Density fitting

.. module:: df
   :synopsis: Density fitting and RI approximation
.. sectionauthor:: Qiming Sun <[email protected]>.

Introduction

The :mod:`df` module provides the fundamental functions to handle the 3-index tensors required by the density fitting (DF) method or the resolution of identity (RI) approximation. Specifically, it includes the functions to compute the 3-center 2-electron AO integrals, the DF/RI 3-index tensor in the form of Cholesky decomposed integral tensor ((ij|kl)=V_{ij,x}V_{kl,x}), the AO to MO integral transformation of the 3-index tensor, as well as the functions to generate the density fitting basis.

The :func:`density_fit` method can utilize the DF method at SCF and MCSCF level:

from pyscf import gto, scf, mcscf
mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='def2-tzvp')
mf = scf.RHF(mol).density_fit().run()
mc = mcscf.CASSCF(mf, 8, 10).density_fit().run()

Once the DF method is enabled at the SCF level, all the post-SCF methods will automatically enable the DF method, for example:

from pyscf import gto, dft, tddft
mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='def2-tzvp')
mf = dft.RKS(mol).density_fit().run()
td = tddft.TDA(mf).run()
print(td.e)

In PySCF, DF is implemented at different level of implementations for different methods. They are summarized in the following table

Methods Fake-ERI Native DF Properties with DF
HF/DFT Yes Yes  
Generlized HF/DFT Yes No  
Relativistic HF Yes Yes  
TDDFT Yes Yes  
RCCSD Yes Yes  
UCCSD Yes No  
RCCSD(T) Yes No  
EOM-CCSD Yes No  
RMP2 Yes No  
UMP2 Yes No  
PBC HF/DFT Yes Yes  
PBC TDDFT Yes Yes  
PBC Gamma-point CCSD Yes Yes  
PBC k-points CCSD Yes No  

Fake-ERI means to mimic the 4-center 2-electron repulsion integrals (ERI) by precontracting the DF 3-index tensor. This is the simplest way to enable DF integrals, although the fake-ERI mechanism may require huge amount of memory also may be slow in performance. It provides the most convenient way to embed the DF integrals in the existing code, thus it is supported by almost every method in PySCF. It is particularly important in the periodic code. Using the fake-ERIs allows us to call all quantum chemistry methods developed at molecular level in the \Gamma-point calculations without modifying any existing molecular code. See also the :ref:`pbc_df` module.

Some methods have native DF implementation. This means the performance of the DF technique has been considered in the code. In these methods, DF approximation generally runs faster than the regular scheme without integral approximation and also consumes less memory or disk space.

When density fitting is enabled in a method, a :attr:`with_df` object will be generated and attached to the method object. :attr:`with_df` is the object to hold all DF relevant quantiles, such as the DF basis, the file to save the 3-index tensor, the amount of memory to use etc. You can modify the attributes of :attr:`with_df` to get more control over the DF methods. In the SCF and MCSCF methods, setting :attr:`with_df` to None will switch off the DF approximation. In the periodic code, all two-electron integrals are evaluated by DF approximations. There are four different types of DF schemes (:class:`FFTDF`, :class:`AFTDF`, :class:`GDF`, :class:`MDF` see :ref:`pbc_df`), available in the periodic code. By assigning different DF object to :attr:`with_df`, different DF schemes can be applied in the PBC calculations.

DF auxiliary basis

The default auxiliary basis set are generated by function :py:func:`pyscf.df.addons.make_basis` based on the orbital basis specified in the calculation according to the rules defined in :data:`pyscf.df.addons.DEFAULT_AUXBASIS`. Specifically, the jkfit basis in the first column is used for Hartree-Fock or DFT methods, and the ri basis in the second column is used for correlation calculations. These optimized auxiliary basis sets are obtained from http://www.psicode.org/psi4manual/master/basissets_byfamily.html If optimized auxiliary basis set was not found for the orbital basis set, even-tempered Gaussian functions are generated automatically.

Specifying auxiliary basis is a common requirement in the real applications. For example, the default auxiliary basis set for the pure DFT calculations may be over complete since it is designed to represent both the Coulomb and HF exchange matrix. Coulomb fitting basis such as Weigend-cfit basis or Ahlrichs-cfit basis are often enough to obtain chemical accuracy. To control the fitting basis in DF method, You can change the value of :attr:`with_df.auxbasis` attribute. The input format of auxiliary fitting basis is exactly the same to the input format of orbital :ref:`gto_basis` set. For example:

from pyscf import gto, dft
mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='def2-tzvp')
mf = dft.RKS(mol)
mf.xc = 'pbe,pbe'
mf.run()  # -109.432313679876
mf = mf.density_fit()
mf.run()  # -109.432329411505
mf.with_df.auxbasis = 'weigend'
mf.run()  # -109.432334646584

More examples for inputing auxiliary basis in the DF calculation can be found in examples/df/01-auxbasis.py.

Even-tempered auxiliary Gaussian basis

The even-tempered auxiliary Gaussian basis is generated by function :func:`aug_etb`:

from pyscf import gto, df
mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='ccpvdz')
print(mol.nao_nr())                            # 28
auxbasis = df.aug_etb(mol)
print(df.make_auxmol(mol, auxbasis).nao_nr())  # 200
auxbasis = df.aug_etb(mol, beta=1.6)
print(df.make_auxmol(mol, auxbasis).nao_nr())  # 338

Here the :func:`make_auxmol` function converts the auxbasis to a :class:`Mole` object which can be used to evaluate the analytical integrals the same way as the regular :class:`Mole` object. The formula to generate the exponents \zeta of the even-tempered auxiliary basis are

\varphi &= r^l \exp(-\zeta_{il} r^2), \quad i = 0..n \\
\zeta_{il} &= \alpha \times \beta^i
:label: etb

The default value of \beta is 2.3. \alpha and the number of auxiliary basis n is determined based on the orbital basis. Given the orbital basis

\chi = r^l \exp(-\alpha_l r^2)

the orbital pair on the same center produces a new one-center basis

\chi \chi' = r^{l+l'} \exp(-(\alpha_l+\alpha_l') r^2)
= r^L \exp(-\alpha_L r^2)

The minimal \alpha_L in all orbital pairs is assigned to \alpha in :eq:`etb`. Then n is estimated to make the largest auxiliary exponent \zeta as close as possible to the maximum \alpha_L. The size of generated even-tempered Gaussian basis is typically 5 - 10 times of the size of the orbital basis, or 2 - 3 times more than the optimized auxiliary basis. (Note the accuracy of this even-tempered auxiliary basis is not fully benchmarked. The error is close to the optimized auxiliary basis in our tests.)

Saving/Loading DF integral tensor

Although it is not expensive to compute DF integral tensor in the molecular calculation, saving/loading the 3-index tensor is still useful since it provides an alternative way, different to the attribute :attr:`_eri` of mean-field object (see :ref:`customize_h`), to customize the Hamiltonian.

In the DF-SCF method, the 3-index tensor is held in the :attr:`with_df` object. The :attr:`with_df` object (see :class:`pyscf.df.df.DF` class) provides two attributes :attr:`_cderi_to_save` and :attr:`_cderi` to access the DF 3-index integrals.

If a DF integral tensor is assigned to :attr:`_cderi`, the integrals will be used in the DF calculation. The DF integral tensor can be either a numpy array or an HDF5 file on disk. When the DF integrals are provided in the HDF5 file, the integral needs to be stored under the dataset 'j3c':

import numpy
import h5py
from pyscf import gto, scf, df
mol = gto.M(atom='H 0 0 0; H 1 0 1; H 0 1 1; H 1 1 0', basis='sto3g')
nao = mol.nao_nr()
with h5py.File('df_ints.h5', 'w') as f:
    f['j3c'] = numpy.random.random((10,nao*(nao+1)//2))
mf = scf.RHF(mol).density_fit()
mf.with_df._cderi = 'df_ints.h5'
mf.kernel()

As shown in the above example, the integral tensor V_{x,ij} provided in :attr:`_cderi` should be a 2D array in C (row-major) convention. Its first index corresponds to the auxiliary basis and the second combined index ij is the orbital pair index. When load DF integrals, we assumed hermitian symmetry between the two orbital index, ie only the elements i\geq j are left in the DF integral tensor. Thus the DF integral tensor should be a 2D array, with shape (M,N*(N+1)/2), where M is the number of auxiliary functions, N is the number of orbitals.

If :attr:`_cderi` is not specified, the DF integral tensor will be generated during the calculation and stored to the file that the attribute :attr:`_cderi_to_save` points to. By default, it is a random file and the random file will be deleted if the calculation finishes successfully. You can find the filename in the output log (when with.verbose > 3, for example:

******** <class 'pyscf.df.df.DF'> flags ********
auxbasis = None
max_memory = 20000
_cderi_to_save = /scratch/tmp6rGrSD

If the calculation is terminated problematically with error or any other reasons, you can reuse the DF integrals in the next calculation by assigning the integral file to :attr:`_cderi`. Overwriting :attr:`_cderi_to_save` with a filename will make the program save the DF integrals in the given filename regardless whether the calculation is succeed or failed. See also the example pyscf/examples/df/10-access_df_integrals.py.

Precomputing the DF integral tensor

The DF integral tensor can be computed without initialization the :attr:`with_df` object. Functions :func:`cholesky_eri` defined in :mod:`df.incore` and :mod:`df.outcore` can generate DF integral tensor in memory or in a HDF5 file:

from pyscf import gto, df
mol = gto.M(atom='N 0 0 0; N 1 1 1', basis='ccpvdz')
cderi = df.incore.cholesky_eri(mol, auxbasis='ccpvdz-jkfit')
df.outcore.cholesky_eri(mol, 'df_ints.h5', auxbasis='ccpvdz-jkfit')

These cderi integrals has the same data structure as the one generated in :attr:`with_df` object. They can be directly used in the DF type calculations:

from pyscf import scf
mf = scf.RHF(mol).density_fit()
mf.with_df._cderi = cderi
mf.kernel()

mf.with_df._cderi = 'df_ints.h5'
mf.kernel()

Approximating orbital hessian in SCF and MCSCF

Orbital hessian is required by the second order SCF solver or MCSCF solver. In many systems, approximating the orbital hessian has negligible effects to the convergence and the solutions of the SCF or MCSCF orbital optimization procedure. Using DF method to approximate the orbital hessian can improve the overall performance. For example, the following code enables the DF approximation to the orbital hessian in SCF calculation:

from pyscf import gto, scf
mol = gto.M(atom='N 0 0 0; O 0 0 1.5', spin=1, basis='ccpvdz')
mf = scf.RHF(mol).newton().density_fit().run(verbose=4)  # converged SCF energy = -129.0896469563
mf = scf.RHF(mol).run(verbose=4)                         # converged SCF energy = -129.0896469563

The approximation to orbital hessian does not change the SCF result. In the above example, it produces the same solution to the regular SCF result. Similarly, when the DF approximation is used with CASSCF orbital hessian, the CASSCF result should not change. Continuing the above example, we can use the :func:`mcscf.approx_hessian` function to change the orbital hessian of the given CASSCF object:

from pyscf import mcscf
mc = mcscf.approx_hessian(mcscf.CASSCF(mf, 8, 11)).run()  # -129.283077136
mc = mcscf.CASSCF(mf, 8, 11).run()                        # -129.283077136

Note

In the second order SCF solver, the order to apply the density_fit and newton methods affects the character of the resultant SCF object. For example, the statement mf = scf.RHF(mol).density_fit().newton() first produces a DFHF object then enable the second order Newton solver for the DFHF object. The resultant SCF object is a DFHF object. See more examples in examples/scf/23-decorate_scf.py

Examples

Program reference

DF class

.. autoclass:: pyscf.df.df.DF
   :members:

df.incore

.. automodule:: pyscf.df.incore
   :members:

df.outcore

.. automodule:: pyscf.df.outcore
   :members:

df.addons

.. automodule:: pyscf.df.addons
   :members:

df.df

.. automodule:: pyscf.df.df
   :members:


df.df_jk

.. automodule:: pyscf.df.df_jk
   :members:

df.r_incore

.. automodule:: pyscf.df.r_incore
   :members:

df.grad.rhf

.. automodule:: pyscf.df.grad.rhf
   :members:

df.grad.uhf

.. automodule:: pyscf.df.grad.uhf
   :members:

df.grad.rks

.. automodule:: pyscf.df.grad.rks
   :members:

df.grad.uks

.. automodule:: pyscf.df.grad.uks
   :members:

df.hessian.rhf

.. automodule:: pyscf.df.hessian.rhf
   :members:

df.hessian.uhf

.. automodule:: pyscf.df.hessian.uhf
   :members:

df.hessian.rks

.. automodule:: pyscf.df.hessian.rks
   :members:

df.hessian.uks

.. automodule:: pyscf.df.hessian.uks
   :members: