hf.cpp line 199

[aromanro]

Hi,

I am implementing a Hartree-Fock program for my blog. I pulled your program to compare results for electron-electron integrals (I am using a different method of calculating them). Of course I still have issues in the implementation, but anyway, I was looking over your code and noticed:

Fp = trimatprod(Xp,F,X);

but then on line 199:

Eigsym es(F,Cc,molorben);

It looks to me as what you intended to have is:

Eigsym es(Fp,Cc,molorben);

I don't know exactly what Eigsym does, I didn't look into it. I suppose it only solves the ordinary eigenvalue problem (that is, not the generalized one), since you don't pass to it the overlap matrix. I think you should pass the transformed Fock operator to it.

I changed the code and executed it for H2O (not your input file, but a configuration I got from a Mathematica Journal article, which is a little bit different) and the results seem to come out a little better.

Thanks,
Adrian

[aromanro]

It seems that now I have the recurrence relations working, I get very close values to the ones hfcxx gets for electron-electron integrals.

Now I get a very close value for the H2O energy with my program compared with hfcxx with the above change.

[ifilot]

Hi Adrian,

Thanks a lot for your comments (also with regard to issue #2). I have not found the time to look into this in detail, but I believe you have indeed found something that is clearly wrong. I'll let you know.

[ifilot]

I have compared the results with PyQuante, with the line

Eigsym es(F,Cc,molorben);

I get the following output for the orbitals of H2:

|**********|**********|**********|**********|
|1         |     -0.9594227156 HT|  Occupied|
|**********|**********|**********|**********|
|[1] H1  (1s)         |              -0.5489|
|[2] H2  (1s)         |              -0.5489|
|----------|----------|----------|----------|

|**********|**********|**********|**********|
|2         |      0.2283480273 HT|   Virtual|
|**********|**********|**********|**********|
|[1] H1  (1s)         |              -1.2115|
|[2] H2  (1s)         |               1.2115|
|----------|----------|----------|----------|

Changing the line to

Eigsym es(Fp,Cc,molorben);

gives

|**********|**********|**********|**********|
|1         |     -0.5782029756 HT|  Occupied|
|**********|**********|**********|**********|
|[1] H1  (1s)         |              -0.5489|
|[2] H2  (1s)         |              -0.5489|
|----------|----------|----------|----------|

|**********|**********|**********|**********|
|2         |      0.6702677791 HT|   Virtual|
|**********|**********|**********|**********|
|[1] H1  (1s)         |              -1.2115|
|[2] H2  (1s)         |               1.2115|
|----------|----------|----------|----------|

I have compared these numbers with the result of PyQuante which confirms that the latter is correct and the former is wrong, in line with this issue.