I would like to plot the XC potential in real space (on a grid) but I am unsure on how to proceed and if there is already an option in Psi4 to do so ? I didn’t find it, but I’m wondering if it can be accessible in a cube file or similarly to plotting the density in real space (see this feature-request, but I don’t think the feature has been released yet).
I’ve also looked at tutorials 4a and 4b but I’m not sure if they answer this question or not.

Note that I actually would like to plot my own V_XC, hence it will be different than the one obtained by doing:

e, wfn = psi4.energy("SVWN", return_wfn=True)
Da = wfn.Da_subset("AO")
Vxc = Da.clone()
Vpotential = wfn.V_potential()
Vpotential.set_D([Da])
Vpotential.compute_V([Vxc])

For the moment, I’m studying the hydrogen molecule in a small basis like 6-31G or cc-pVDZ. A 1D plot along the bonding axis is what I would like to do.

This depends on what you mean by the XC potential. For a LDA type functional, you can compute the value from the density, which Psi4 can already dump in a cube file.

I think that what compute_V computes is actually a Fock matrix contribution: Vpotential is a matrix with elements <u|F|v>. One would have to think if it is possible to project this into real space…

Yes, for a LDA-type functional I could just use the density directly.
But I am indeed interested about projecting any V_xc matrix (the Fock matrix contribution) into the real space. I believe this should be possible, because the matrix is expressed in a given basis (in my case, I use orthogonal atomic orbitals).

So I would guess the following steps:

compute my new V_xc in the OAO basis

Transform back this matrix into the AO basis

As the AO basis functions are GTOs, I think I can define V_xc on a grid, similarly as the density matrix in the AO basis as also projected into the real space.

Step 1 and 2 are done but step 3 is much more complex, at least to me. Is my reasoning wrong and the V_xc matrix cannot be projected in real-space for some specific reason ?

Thanks for the link. However, I believe I want something different. As far as I know, inverse KS theory is able to find a V_XC potential from a given density. But I already have a potential (expressed in the AO basis), I just want to know how if it is possible to plot it in real space, and if yes how can I do that ?

Right, sorry; I had forgotten to reply and apparently forgotten what the issue was…

It might be as simple as multiplying F with the matrix of basis function values X on the grid: V = XFX^T. However, there is no guarantee that V will be local; in fact, exact exchange is by definition non-local, V is therefore an Ngrid^2 quantity.

In addition, often you have to introduce factors of S^{-1} so that projecting back on to the basis set gives you the correct result. In fact, this is exactly what it would appear: if you have V given on the grid as above, and you want to project it back into the basis set, you would introduce another pair of Xs on the left and right, as well as do quadratures over the two grids over a weight function. You will get two terms XX^T which give you S integrated over the grid. Therefore, it appears that the correct definition is actually

Thanks ! That’s indeed exactly what I want. I agree that V will be non-local in general, that’s fine. I assume S^{-1} is the inverse of the overlap matrix in the AO basis, but do you know how I can extract that matrix X in psi4 ? I believe it is exactly the “PHI” in this tutorial?