Overwrite density matrix with new basis set(s)


I would like to calculate the electrostatic potential and field on a grid from my own, custom density matrix. This is very similar to Cube file from arbitrary density matrix, but with one big difference. I want to specify new basis sets for my density matrix. I am representing my density in an auxiliary basis, rho = sum(c_i * phi_i), so my density matrix will not be square. It will have dimension (num_aux_ao, 1), using the zero_ao_basis_set.

Is there any way psi4 can do this right now? It looks to me like the cubeprop code or the analogous code in oeprop would be a decent place to start.


How is that a density matrix? It can be an electron density, but what you are describing does not at all sound like an density matrix.

@jmisiewicz, it’s a density, and it’s a matrix, so that makes it a “density matrix”, right?


Kidding aside, I think there must be a way to get psi4 to give me an electrostatic potential from my density. Like I say above, I think this looks like defining a density “matrix” that is really just a vector.

There are two pieces in psi4 that make me think this might be possible.

  1. In the cubeprop code, (CubicScalarGrid::add_esp in particular) psi4 is already doing density fitting before computing the electrostatic field on a grid. Since what I have is directly analogous to the auxiliary basis set density here, I think there could be a way to make this code do what I want.
  2. In the oeprop code, I’m not sure that it is requiring a square matrix as input. In ElectrostaticInt::compute_pair it looks to me like it is taking two independent basis sets. I could be wrong about this. Let me know if I am.

If it is not currently possible with #1 or #2, I am willing to do some hacking to get my use case to work. I really would like to avoid writing my own integral code here…


libmints/oeprop.cc uses the 1RDM, which is not what you have, so I suspect libcubeprop is going to be of more use to you. I’m not familiar with this library, so I can’t be of further assistance.

Sounds good. Thanks, @jmisiewicz.

@andysim, @dgasmith or others who wrote the cubeprop library, any ideas on how I might be able to do this? I’m happy to get my hands dirty, assuming some code modifications are necessary, just looking for some general guidance.

Another idea I have is to project my density in the auxiliary basis back into the original orbital basis. It’s not obvious to me that density-fitting is invertible, though…

Sorry for the spam, but yet another idea here.

I think ao_potential can be used with a mixed basis.

Overloaded function.

ao_potential(self: psi4.core.MintsHelper) -> psi4.core.Matrix
AO potential integrals

ao_potential(self: psi4.core.MintsHelper, arg0: psi::BasisSet, arg1: psi::BasisSet) -> psi4.core.Matrix
AO mixed basis potential integrals

If I interpret this right. I think I can call something like this:

potential_mat = mints.ao_potential(aux_basis, zero_basis)

Does this look right? If so, I just need to figure out how to calculate esp on a grid from potential_mat.

Post 1:

For technical help on this level, you’re better off trying the Slack channel.

Post 2:

I don’t understand how that would help. Electron densities are fundamentally not convertible back to 1RDMs, and you need to use grid-based techniques. The functional form of your electron density is beside the point.

Post 3:

No, that is most definitely wrong. That function has absolutely nothing to do with an electrostatic potential. Those are electron-nuclear attraction integrals between two orbitals. If you tried to pass in a zero element basis set, there would be zero columns in the resulting matrix.

This is just density fitting. You expand the density matrix Puv in terms of the auxiliary basis set as c_B = \sum_{uv} Puv (uv|A) (A|B)^{-1}; this is then the expansion of the electron density. You get the electrostatic potential as a linear combination of the potentials of the auxiliary functions (easily evaluated analytically), and the field is the gradient of the potential integrals.

Hi, @susilehtola!

This is exactly right. I am doing density fitting.

You get the electrostatic potential as a linear combination of the potentials of the auxiliary functions

Precisely. I just want to know if psi4 already has this coded up anywhere. Chalk it up to a combination of laziness in not wanting to write these integrals up myself and worries about staying consistent with how psi4 defines its basis sets.