.. _user_df: Density fitting (DF) ******************** *Modules*: :py:mod:`pyscf.df`, :py:mod:`pyscf.pbc.df` .. module:: df :synopsis: Density fitting and RI approximation .. sectionauthor:: Qiming Sun . Introduction ============ Density fitting (DF), sometimes also called the resolution of identity (RI) approximation, is a method to approximate the four-index electron repulsion integrals (ERIs) by two- and three-index tensors. :cite:`Dunlap00PCCP` In DF, the atomic orbital (AO) product space is expanded in terms of an auxiliary basis set. PySCF has extensive support for using DF in both molecular and crystalline calculations. This section covers the basic usage of DF in molecular applications. See :ref:`user_pbc_df` for the use of DF in periodic calculations. Using DF in SCF calculations ============================ DF is not used by default when initializing a :class:`SCF` or :class:`MCSCF` object, but can be invoked via the :func:`density_fit` method as follows:: 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() # output: -108.943773290737 mc = mcscf.CASSCF(mf, 8, 10).density_fit().run() # output: -109.077600181849 Alternatively, one can initialize a :class:`DF` object and overwrite the :attr:`with_df` attribute of the :class:`SCF` object:: from pyscf import gto, scf, df mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='def2-tzvp') mydf = df.DF(mol).build() mf = scf.RHF(mol).density_fit() mf.with_df = mydf mf.kernel() # output: -108.943773290737 Setting :attr:`with_df` to ``None`` turns DF off. Using DF in post-SCF calculations ================================= A post-SCF calculation will use DF automatically if the preceding SCF calculation is performed using DF. The example below performs a TDDFT/TDA calculation using DF:: 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, xc="b3lyp").density_fit().run() td = tddft.TDA(mf).run() print(td.e) # output: [0.27294247 0.28837975 0.2883798 ] Alternatively, one can use the :func:`density_fit` method provided by most post-SCF classes in the same way as in a SCF calculation described above. This is useful when one wants to use DF to accelerate only the post-SCF calculations or to use a different auxiliary basis for the post-SCF calculations than that used by SCF. The example below performs a DF-MP2 calculation using the ``cc-pvtz-ri`` basis :cite:`Weigend02JCP` on top of a SCF calculation that does not use DF:: from pyscf import gto, scf, mp mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='cc-pvtz') mf = scf.RHF(mol).run() mmp = mp.MP2(mf).density_fit(auxbasis="cc-pvtz-ri") e_corr = mmp.kernel()[0] print(e_corr) # output: -0.4288734824009443 .. _choice_of_auxbasis: Choice of auxiliary basis ========================= General considerations ---------------------- The choice of auxiliary basis depends on the AO basis. By default, PySCF uses a pre-defined auxiliary basis set optimized for the used AO basis set, if one exists. This includes many commonly used AO basis sets in electronic structure calculations, e.g., the Ahlrichs `def2` family, :cite:`Hellweg07TCA,Weigend98CPL` and the Dunning `cc` family, :cite:`Weigend02PCCP`. When a pre-defined auxiliary basis set is not available, an even-tempered basis (ETB) set is generated by the following rule .. math:: \varphi_i = r^l \exp(-\alpha_l \beta^i r^2), \quad i = 0, 1, \cdots, n_l where both :math:`\alpha_l` and :math:`n_l` are determined automatically from the AO basis. Specifically, :math:`\alpha_l` is set to the *smallest* exponent of the AO products with angular momentum :math:`l`, and :math:`n_l` is chosen so that the *largest* exponent of the AO products with angular momentum :math:`l` lies in between :math:`\alpha_l \beta^{n_l-1}` and :math:`\alpha_l \beta^{n_l}`. The parameter :math:`\beta` controls the "density" of the fitting functions: smaller :math:`\beta` gives an ETB of larger size and *vice versa*. The default is :math:`\beta = 2.0`. The user can overwrite the default choice of the auxiliary basis by setting :attr:`auxbasis` upon initialization or at a later stage:: from pyscf import gto, scf, df mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='cc-pvdz') mf = scf.RHF(mol).density_fit(auxbasis="weigend") mf.kernel() # -108.910953335055 mf.with_df.auxbasis = "cc-pvdz-jkfit" # this is the default mf.kernel() # -108.913710743723 mf.with_df.auxbasis = df.aug_etb(mol, beta=1.7) # ETB with beta = 1.7 mf.kernel() # -108.914059329528 The user can print the used auxiliary basis set at any time by:: print(mf.with_df.auxmol.basis) More examples on using customized auxiliary basis can be found in :source:`examples/df/01-auxbasis.py`. Special considerations for DFT ------------------------------ For DFT calculations with pure exchange-correlation functionals (i.e., LDA and GGA), the default auxiliary basis, which is designed for fitting both the Coulomb and the Hartree-Fock exchange integrals, may be unnecessarily large. We recommend using the ``def2-universal-jfit`` basis :cite:`Weigend06PCCP` for a more cost-effective choice as shown in the following example:: from pyscf import gto, dft mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='def2-tzvpp') mf = dft.RKS(mol, xc="pbe").density_fit().run() # -109.432329411505 print(mf.with_df.auxmol.basis) # default: def2-tzvpp-jkfit print(mf.with_df.auxmol.nao_nr()) # 154 aux basis functions mf.with_df.auxbasis = "def2-universal-jfit" # same as "weigend" mf.kernel() # -109.432334646585 print(mf.with_df.auxmol.nao_nr()) # 98 aux basis functions .. _save_reuse_df_tensors: Saving and reusing DF tensors ============================= A key result of the :class:`DF` class is the Cholesky decomposition of the ERIs (CDERIs) in terms of the 3-index tensor :math:`d_{L\mu\nu}`, .. math:: (\mu\nu|\lambda\sigma) = \sum_{L} d^*_{L\mu\nu} d_{L\lambda\sigma} The :func:`build` method of :class:`DF` computes the CDERIs. By default, the 3-index tensor defining the CDERIs is discarded once the calculation finishes. Sometimes it is useful to save it to disk for re-use in later calculations. This can be achieved by specifying a HDF5 file by setting :attr:`_cderi_to_save` either at the SCF stage:: mf = scf.RHF(mol).density_fit() mf.with_df._cderi_to_save = 'saved_cderi.h5' mf.kernel() or by initializing a :class:`DF` object separately:: mydf = df.DF(mol) mydf.auxbasis = df.make_auxbasis(mol) mydf._cderi_to_save = 'saved_cderi.h5' mydf.build() The saved DF tensor can be used later by setting :attr:`_cderi` to the HDF5 file:: mf = scf.RHF(mol).density_fit() mf.with_df._cderi = 'saved_cderi.h5' mf.kernel() More examples on saving and using DF integrals and tensors can be found in :source:`examples/df/10-access_df_integrals.py`, :source:`examples/df/11-access_df_tensor.py`, and :source:`examples/df/40-precompute_df_integrals.py`. Advanced examples ================= More advanced examples of using the :mod:`df` module include * Computing the :math:`\mathbf{J}`-matrix of DFT in an I/O-free manner: :source:`examples/df/11-get_j_io_free.py`. * Using DF integrals to define the Hamiltonian of a CASSCF calculation: :source:`examples/df/40-custom_df_hamiltonian.py`. * Generating analytical gradients for DF integrals: :source:`examples/df/41-df_integrals_gradients.py`. * Customizing the :func:`get_jk` method of a SCF class using DF: :source:`examples/df/42-overwrite_get_jk.py`.