Team,

I am working on a tool that requires as it's input an electron density on a grid. This was easy to achieve for the neutral and ionic states for given elements by using the numint.eval_rho function. However, when I tried to move to excited states, calculated using TDDFT, the densities being produced using numint.eval_rho2 function varied significantly between runs, and gave results that were non-sensical.

How should I calculate the electron density on a grid for an excited state of an element?

Is the numint.eval_rho2 function the right one to use with a TDDFT function in order to extract the electron density? Here is the code I am using. I previously tried the other eval_rho functions without any success. Please let me know if I am on the wrong path here.

Honestly, I am hoping that there is a tool that I am missing or an approach that I have not found yet.

Thank you in advance!

mf = dft.UKS(mol)
mf.xc = 'LDA_X,LDA_C_PZ'
mf.kernel()  

fin = pyscf.tdscf.uks.TDA(mf)
fin.nstates = 1
fin.kernel()

grids = dft.gen_grid.Grids(mol)
grids.level = 2 
grids.build()

coords = grids.coords

weights = grids.weights

ao_value = numint.eval_ao(mol, coords, deriv=1)

rho1, dx_rho1, dy_rho1, dz_rho1 = numint.eval_rho2(mol, ao_value, fin._scf.mo_coeff[0,:,:], fin._scf.mo_occ[0,:], non0tab=None, xctype='GGA')

rho2, dx_rho2, dy_rho2, dz_rho2 = numint.eval_rho2(mol, ao_value, fin._scf.mo_coeff[1,:,:], fin._scf.mo_occ[1,:], non0tab=None, xctype='GGA')