Hi, I’m new to the forum.
I have a question external to psi4, but I hope someone can help me. I’m programming a small HF code which handles 2-el integrals up to [pppp] shells for now. I managed to reproduce exactly reported integrals for STO-3G methane, and also the SCF energy. However, when I changed to an acetilene molecule (with p functions on more than one center), I got a deviation on the energy. Also If I do 6-31g calculations on CH4 I get a value which is off (~-46 Ha versus ~-40 Ha reference value), and if I do 6-31g on acetile I get around -9000 Ha.
For these types of calculations I can only compare the energy or the one-electron matrices (which are correct to about 1e-12). I’m almost sure the problem lies in integrals which contain [pp **] shells with p functions which are different (either different centers or different exponents), and my bug vanishes for highly symmetric cases, but after weeks I cannot find the problem.
My question is if someone could provide me with a list of integrals for some small molecule containing more than one C, N or O atom with an STO-3G or 6-31g basis. (I could also code some other basis, but for now I can only handle s and p type functions)
I know it can be done with psi4. I have downloaded it but still wasn’t able to use it because I’m having problems installing a new version of python. I don’t know hoy much time it will take me to solve the python issue and learn to use psi4, and this integral bug is breaking my head.
If someone could help me I would deeply appreciate it.
Best regards,
Nicolas
You might want to check out the Psi4Numpy repo here: https://github.com/psi4/psi4numpy
I don’t have a link on my phone but Psi4 is conda installable which gets a complete binary package. This should be in the installation part of our manual.
I would recommend debugging against LIBCINT. it is easy to install it, however a bit tedious to operate if you are not a pro at C (PySCF offers a more friendly interface to the library but I never managed to get around it).
That way you can check values for individual <ab|cd> integrals with arbitrary GTOs (be mindful of normalisation and L-ordering). Should be faster to find the bug if it is indeed in the ERI algorithm. The full J/K matrixes also involve a sum over occupied states of appropriate spin… the bug may be there too.
In the meantime, here are some data for some small organics:
There are 6 molecules from the Rupp/Ramakrishnan dataset. Geometries are in the input files.
The .dat files are the overlap, Fock (Fa) and one-electron (H) matrixes.All matrix elements below 1.0e-16 were zeroed for sanity.
Calculations are direct RHF SCF with STO-3G basis from PSI4. I am not sure what exactly is in the Fock matrix… I printed out wavefunction.Fa(), so I guess it is H+J-K for the spin-up electrons only. Devs shall correct me if wrong. You should be able to get the (J-K) matrix but probably not the individual parts. I am also adding now the explicit ERIs obtained from the MintsHelper class in PSI4. Let me know if you want more molecules. Hope this helps.
Just out of curiosity… how are you dealing with the ERI? Are you using some well established scheme or some new method?
Felix, thank you very much for the answer and data. It’s really useful.
When I posted the other day I had a large difference in energy. Finally I found out that the problem was not in the ERIs but on the way I was using the 8-fold permutational symmetry of the ERIs to do only the unique integrals. I constructed the ERI(a,b,c,d) matrix using this symmetry, and when I switched back to doing the full matrix, there were about 4 elements (in 17^4) that had some deviations. I still don’t know why is that, but started doing all the integrals, and then my energy was almost correct. I was about to post an answer here, but noticed that instead of -40.18… Ha for 6-31g methane obtained with Orca, I got -40.20… Ha. This small difference I still can’t figure out. So there is still some problem.
For the ERIs I’m using Obara-Saika derived methods, something similar to Head-Gordon Pople. I use the vertical recurrence to obtain the [L s s s] integrals (up to [g s s s] integrals with one algorithm and up to [L = 10 s s s] with another algorithm, slightly slower but allowing arbitrary angular momentum). I use the electron transfer recursion to obtain [La s Lc s] integrals, then I contract and then use horizontal recursion to obtain the full integral matrices for each shell. The code is vectorized, and works on all cartesian components simultaneously.
I do permutations to get all shells into the order I just described, and I think the permutations are correct.
Once again, thank you very much,
I’ll post when/if I solve the problem.
Best,
Nicolas