1#!/usr/bin/env python 2# 3# Author: Qiming Sun <osirpt.sun@gmail.com> 4# 5 6from pyscf import gto 7from pyscf import scf 8 9''' 10Density fitting method by decorating the scf object with scf.density_fit function. 11 12There is no flag to control the program to do density fitting for 2-electron 13integration. The way to call density fitting is to decorate the existed scf 14object with scf.density_fit function. 15 16NOTE scf.density_fit function generates a new object, which works exactly the 17same way as the regular scf method. The density fitting scf object is an 18independent object to the regular scf object which is to be decorated. By 19doing so, density fitting can be applied anytime, anywhere in your script 20without affecting the exsited scf object. 21 22See also: 23examples/df/00-with_df.py 24examples/df/01-auxbasis.py 25''' 26 27mol = gto.Mole() 28mol.build( 29 verbose = 0, 30 atom = '''8 0 0. 0 31 1 0 -0.757 0.587 32 1 0 0.757 0.587''', 33 basis = 'ccpvdz', 34) 35 36mf = scf.density_fit(scf.RHF(mol)) 37energy = mf.kernel() 38print('E = %.12f, ref = -76.026744737355' % energy) 39 40# 41# Stream style: calling .density_fit method to return a DF-SCF object. 42# 43mf = scf.RHF(mol).density_fit() 44energy = mf.kernel() 45print('E = %.12f, ref = -76.026744737355' % energy) 46 47# 48# By default optimal auxiliary basis (if possible) or even-tempered gaussian 49# functions are used fitting basis. You can assign with_df.auxbasis to change 50# the change the fitting basis. 51# 52mol.spin = 1 53mol.charge = 1 54mol.build(0, 0) 55mf = scf.UKS(mol).density_fit() 56mf.with_df.auxbasis = 'cc-pvdz-jkfit' 57energy = mf.kernel() 58print('E = %.12f, ref = -75.390366559552' % energy) 59 60# Switch off density fitting 61mf.with_df = False 62energy = mf.kernel() 63print('E = %.12f, ref = %.12f' % (energy, scf.UKS(mol).kernel())) 64