from pyscf import gto,scf,dft

# 1. Define System

mol         = gto.M(atom   = 'C 0 0 0; O 0 0 1.5',

                    basis  = 'ccpvdz',

                    charge = 0,

                    spin   = 0)

# 2. Do Hartree-Fock Calculation & Obtain HF Density

myHF        = scf.RHF(mol)

myHF.kernel()

dm_HF       = myHF.make_rdm1()

# 3. Put HF Density into DFT Hamiltonian

myDFT       = dft.RKS(mol)

myDFT.xc    = 'pbe'

eHFDFT      = myDFT.energy_tot(dm_HF)