Plotting XC potential in real space


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 ="SVWN", return_wfn=True)
Da = wfn.Da_subset("AO")
Vxc = Da.clone()
Vpotential = wfn.V_potential()

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:

  1. compute my new V_xc in the OAO basis
  2. Transform back this matrix into the AO basis
  3. 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 ?