Obtain the energy from DFT at arbitrary density matrix

[TimGould]

Hi All,

In Psi4/Numpy it is possible to easily compute a potential (I assume xc remainder potential) for an arbitrary density matrix using VBase, per the example here.

Is there a similarly easy way to calculate the xc energy?

EDIT: If there is a similar way to obtain the total energy of a HF, DFT or hybrid-DFT approximation at arbitrary density matrix that would be as good or better!

Obviously I can do this using blocks, a grid, and direct calls to XC routines. But the internal C++ code will be faster and (hopefully) easier.

All the best,

Tim

[hokru]

Hi Tim,

This would be an interesting feature to have. I think ADF has sth like this. Changing the DFT functional after the SCF.

For the Psi4numpy route the DFT tutorials show the principles how to go from density matrix to xc energy: https://github.com/psi4/psi4numpy/blob/master/Tutorials/04_Density_Functional_Theory/4b_LDA_kernel.ipynb

On the C++ side I only see a function to compute the potential, not only the energy, although it does not seem too difficult to add.

Also tagging @dgasmith, if he has further insights.

[dgasmith]

Similarly you can use the built in C++ routines:

import psi4

# Obtain a SCF Wavefunction
mol = psi4.geometry("He\nsymmetry c1")
scf_e, scf_wfn = psi4.energy("PBE/sto-3g", return_wfn=True)

# Grab a few matrices
D = scf_wfn.Da().clone()
V_out = scf_wfn.Da().clone()
V_out.np[:] = 0

Vpotential = scf_wfn.V_potential()

# Set arbitrary density
D.np[:] = 1

# Compute V which evaluates the energy as well
Vpotential.set_D([D])
Vpotential.compute_V([V_out])
print(V_out.np)

# Actual values exist on this dictionary
print(Vpotential.quadrature_values()["FUNCTIONAL"])

Building out a Vpotential can be done without the wavefunction, but it’s slightly tricky to do so. In principle the signature is here, but (and I am recalling) there is some sort of issue:

• The SCF C++ call
• The export is here

Maybe I’m wrong, that looks like enough to stitch it together from the Python side. I don’t have time to try it out, but let me know if it works!

[TimGould]

Perfect. Can indeed confirm that seems to work for obtaining the energy with varying density matrix (which was my main goal).

In future I’ll see if I can work out a sneaky way to use it for changing the functional itself (e.g. for density-driven / functional-driven analysis and density corrected DFT). But that’s lower priority for me right now.

Thanks for your help,

Tim