pkozlows / fci

0 stars 0 forks source link

How to debug my program by orders of mag? #22

Closed pkozlows closed 1 year ago

pkozlows commented 1 year ago

I am gating a fci e of -7.743700621146326. HF e is -7.739373948970316. Correct fci e is -7.8399080148963369. So, it seems like the energy that I'm getting is not stabilized enough. What did you do when you were writing your program for designating between whether the problem was with the pace factor, spin, or r of summation? I notice that correcting my face factor function did not change the fci energy was getting by much. Did you debug by considering the changes of each of these with respect to how much they changed fci e, i.e. sim to an order of magnitude test, or something else?

Walter-Feng commented 1 year ago

It was more like a trial & error for me when my number is approaching the correct energy but not there yet.

The only way I can think of is to make sure that each section is absolutely correct, from the easiest to the hardest. For example, phase factor will be an easier step to guarantee it is working as expected, via unit test of e.g. random pair of determinants and manually check if it matches your derivation.

Ah shit someone tweeting like a machine gun ...

Then it will boil down to matrix element with one-electron / two-electron difference. Between these, two-electron will be the easiest, as there is only 1 two-electron integral element being involved in the calculation.

After this, I would think if I still understand the one-electron difference case correctly. Was I wrong again about spins? How is it written in CI string? In what case the contribution will be zero? etc. And I was able to identify the problem with spin, corrected, and energy was finally matching.

No I didn't take any drafts during this process, all within my brain. So maybe you can try.

Is it dragging you too much? maybe it's about time to show my solution.

pkozlows commented 1 year ago

I feel like am very close to correct fci e. just for shits and giggles, i commented out a part of my two diff case, and got correct fci e to 2 decimal places, which i haven't gotten before, so i feel that the issue must be there. i understand that for [mp|nq] m and p, and n and q have to have the same spin, but I'm heading trouble wraping my head around fact that spin has to be preserved for the total excitation.

You said this:

Between these, two-electron will be the easiest, as there is only 1 two-electron integral element being involved in the calculation.

idk how to get from condon rule which has two ints -> taking account spin and simplifying down to one. that is I had been able to get this for in my deriviation but Stuck where to go from here regarding the what simplifications/factors are introduced by the delta Functions.

https://github.com/pkozlows/fci/blob/45699e0bafdaca79f67a1f669b066b08443f1d93/latex/two_differences/two_differences.pdf

pkozlows commented 1 year ago

Never mind. I think I have a decent handle on the theory, But I still don't understand this part.

Then it will boil down to matrix element with one-electron / two-electron difference. Between these, two-electron will be the easiest, as there is only 1 two-electron integral element being involved in the calculation.

Especially

there is only 1 two-electron integral element being involved in the calculation

Walter-Feng commented 1 year ago

Ah ok I forgot there can be both the coulomb / exchange terms existing, so actually 2 instead of 1 in some cases

but still there's no summation, so I believe my statement is roughly valid.

Walter-Feng commented 1 year ago

and I understand my choice of words are confusing, so more like

< 0 | H | p -> q, r -> s> = A - B

where A and B are two scalars.

pkozlows commented 1 year ago

My fci e is at -7.8362021822923005 and correct e is -7.8399080148963369, so not there yet, but thinking this means im closer than before? i don't see how i can debug this further rn, but then again maybe i will get some inspiration with time? i want to give it a few more days, and then i will ask for a solution if its ok :)

Walter-Feng commented 1 year ago

Yeah it looks closer, I would guess your phase factor is done correct, and you have treated most of the excitations correct, maybe one or two case where you did not have all the terms or you added something extra. You may want to check two-electron excitation cases more carefully, I would guess there is one term scientifically wrong in one of the cases, considering the magnitude of the error.

pkozlows commented 1 year ago

I think I am treating the two election excitation cases correctly now, but my energy is diverging from almost correct one that I had before. am I thinking about something very incorrectly theoretically or could it be something much different? https://github.com/pkozlows/fci/blob/7c970bf6c484ed6a9d0ccbb68801dd1af4f6da06/latex/two_differences/two_differences.pdf when you said

You may want to check two-electron excitation cases more carefully, I would guess there is one term scientifically wrong in one of the cases, considering the magnitude of the error. and specifically one of the cases I realized that it would be correct to split them into multiple cases, which I wasn't do wing before. but this is diverging my e from the correct value, so I am getting confused.

pkozlows commented 1 year ago

look like i dunno markdown very well yet lol

Walter-Feng commented 1 year ago

Garnet: "I am very confused"

But I would still say I can't make it through with the differentiation of two cases - the grammatical error exactly lies where I believe is important.

I would guess you are describing case 1 as both excited electrons having the same spin. Then it looks reasonable to me.

Case 2 - why do you think the expectation value will kick in? Expectation values are scalars, while you are calculating matrix elements.

pkozlows commented 1 year ago

I am sincerely sorry about overlooking the gramatical errors in the last draft. Next time, I need to proofread my drafts more before I send them to you. I think this should make more sense. https://github.com/pkozlows/fci/blob/30ad289fb3de9b73311466e11bcf749220791fb4/latex/two_differences/two_differences.pdf

Walter-Feng commented 1 year ago

Indistinguishablility is already implemented in second quantization - you don't need to care about it.

You already declare you are talking about "different spins" in case 2, then at least there won't be "the same or different" cases affecting. Further, superposition only means a state being a linear combination of two basis functions, but you are calculating a matrix element between two basis functions, so no superposition.

I also get criticized thousand and one times from Garnet about this, you need proof to validate your statements, it can't be "I think this is right, so this is right". At least, you didn't get correct energy right?

pkozlows commented 1 year ago

https://github.com/pkozlows/fci/blob/2c33f952a1a40635255546468d817291a127b2a1/latex/two_differences/two_differences.pdf

Walter-Feng commented 1 year ago
  1. The creation/annihilation operators can be permutated arbitrarily to calculate the matrix elements - this means permutation does not affect the state or the Hamiltonian. Otherwise all the derivations will be wrong - which is not actually a proof but some evidence. Commutation rules / anti-commutation rules explicitly reflect the indistinguishability.

  2. The quantum mechanics will only be states, observables, and linear transformations. So instead of "electron a with spin upward excited from orbital m to n and electron b with spin ..." blabla, we should be saying "we are calculating matrix element between two states, between which there is only two-electron difference of the occupation, which is m-alpha spin orbital and blabla". This does not break any rule about electrons being indistinguishable, and the property is well preserved in the simplification of second-quantization representation. In harmonic oscillators there is nothing different.

  3. An extreme case, you have helium, whose ground state is 1s up and 1s down. you performed a XPS, you vanished both electrons (i.e. to a galaxy far far away). you performed a (1 alpha + 1 beta) excitation, and I don't think this is something prohibited.

pkozlows commented 1 year ago

giving the two electron differences another bash :) https://github.com/pkozlows/fci/blob/7c416ce3af2d3612b96b71b0a575ee2f91e0724b/latex/two_differences/further_two_differences.pdf

Walter-Feng commented 1 year ago

Then we are on the same page that we can have three cases,

  1. alpha, alpha -> alpha alpha
  2. alpha, beta -> alpha, beta
  3. beta, beta -> beta, beta for the excitation. 1 and 3 are similar cases without any special mathematics, and we both agree that 2 may be something different.

Then 2 you can derive explicitly what it should equal to, without any inference from expectation values or modification to the definition of delta functions. Nothing about superposition, nothing about mean values.

pkozlows commented 1 year ago

https://github.com/pkozlows/fci/blob/bcd46b69c2c624dc676ac51d667a528b707de88f/latex/two_differences/two_differences_v3.pdf

Walter-Feng commented 1 year ago

I slightly forgot what [p] refers to, I vaguely recall it's spin.

Then why would you think these two sub-cases can both stand? Say you have,

m -> alpha n->beta p->alpha q->beta

then you can only have excitation m->p and n->q, which corresponds to case 2.1. You don't alter the spin in m->alpha to q->beta.

pkozlows commented 1 year ago

I guess I'm not just not understanding what is preventing all 4 cases:

  1. m -> alpha n->beta p->alpha q->beta
  2. m -> alpha n->beta p->beta q->alpha
  3. m -> beta n->alpha p->beta q->alpha
  4. m -> beta n->alpha p->alpha q->beta

I slightly forgot what [p] refers to, I vaguely recall it's spin. you were correct about this. You don't alter the spin in m->alpha to q->beta. I'm not completely understanding what you mean by this.

Walter-Feng commented 1 year ago
  1. only [ m p | n q ] is contributing, [ m q | n p ] vanishes
  2. only [ m q | n p ] is contributing, [ m p | n q ] vanishes
  3. only [ m p | n q ] is contributing, [ m q | n p ] vanishes
  4. only [ m q | n p ] is contributing, [ m p | n q ] vanishes

So I don't think there will be two terms.

Walter-Feng commented 1 year ago

I would guess our misconception lies in the determinant basis. I think about FCI in spin orbitals, so it comes natural to me there is only one term contributing in "one alpha + one beta excitation" without treating the 4 sub-cases. Maybe you treated first the orbitals and then the spins, so you have the idea of 4 sub-cases.

This may come slightly weird in phase factors, as phase factors can kick in naturally in spin orbital description, but it might be different when it is done per orbital then per spin. This mainly involves your minus sign in 2.2.

pkozlows commented 1 year ago

https://github.com/pkozlows/fci/blob/7d853fc0ba1436646301d8a387fe38170c222844/latex/two_differences/two_differences_v4.pdf this is the determinant basis generation that I reference in latex file that I think is in term of spin orbitals, but I could be wrong. https://github.com/pkozlows/fci/blob/7d853fc0ba1436646301d8a387fe38170c222844/%20main.py#L22

Walter-Feng commented 1 year ago

It's relatively hard to unravel the misunderstanding we have here, let me briefly describe how I did.

  1. determinant string is unique and should be unique
  2. given 1, a phase factor is unique for arbitrary FCI matrix element
  3. given 1, when m, n -> p, q has mixed spins, with each letter representing spin-orbital, it can be sorted such that m and p always have alpha spins, and n & q always have beta spins.

I would guess you are aligning the spin orbitals as 0a, 0b, 1a, 1b, 2a, 2b, 3a, 3b... with a for alpha spin and b for beta spin, numbers for spatial orbital indices.

say you have 0a, 0b -> 1a, 2b excitation, so you will write [01 | 02] for [mp|nq]. If we set [m]=a, [n]=b, you will perform 2 + 4 + 0 + 0 = 6 permutations, so phase factor = 1. If we set [m]=b, [n]=a, you are exchanging both pairs of operators on the two sides, so 2 + 4 + 0 + 0 + 2 = 8 permutations, so phase factor = 1.

So either way of defining your order of spin when converting matrix element to integral value, the phase factor should not be affecting your result, as phase factor is unique given the definition of determinant string is unique.

I guess here is something important - you should be explicitly defining how you align the alpha/beta spins after you moved the operators involved in the excitation next to Hamiltonian operator. If you define 0a, 0b, 1a, 1b, ...

then bubble sort will align as < beta alpha | H | alpha beta > in mixed spin case. You should not be having a negative sign in any matrix elements.

Note that I am trying to help you identify the bug through our discussion of theory in two-electron integrals. It can be that you are correct, and the problem lies somewhere else; it can be your theory is correct, but your code is not a faithful representation of your theory; it can be your number actually matches the reference value, but some operations are not optimized such that the error accumulated to 10^-3 level.

pkozlows commented 1 year ago

I finally got correct anergy. I feel like a little kid who got a new toy! Thanks for being patient with me :) I'm getting a 404 page not found error for https://sunqm.github.io/pyscf/tutorial.html. do you know the status of how I can learn pyscf?

Walter-Feng commented 1 year ago

Congrats! Would you mind reminding what was the actual problem?

Next step will be writing Davidson diagonalization, and then the Knowles-Handy algorithm that didn't take much time from me.

Do you need to learn pyscf? You can follow the official link to have some quick start.

pkozlows commented 1 year ago

Separating my determinant into alpha and beta parts was definitely helpful. The final thing that I was missing is that I was being boneheaded, and I was implementing the spin exchange part of the no differences incorrectly, which HF assert did not catch because the alpha and beta parts are the same, if that makes sense.

Install PySCF and go through the tutorial at the following URL: "https://sunqm.github.io/pyscf/tutorial.html" to get familiar with the package, so that you can access the integrals for other systems as you want.

I just saw this in the instructions file, and so I was not able to find it. Do I want to be able to access integral for other systems at this point?

I have used PySCF in the past a little bit, but I am not so familiar with it now. is it worth spending sometime now to learn it more?

I just have one final question about my program. In order to produce the FCI matrix, I am looping through every element in the basis twice over, so this is a little bit redundant, as most pairs are getting evaluated twice. After taking a look at what you did in your program, I realize that I could make this more efficient. However, if I am constructing the correct FCI matrix, should this be good enough for now?

pkozlows commented 1 year ago

I guess I should add that the time from the start of my program to getting and FCI anergy is taking just over 2.5 s, which I don't think is too horrible.

Walter-Feng commented 1 year ago

Sounds fantastic - considering that you are doing direct diagonalization plus not using sparse matrix form.

No you don’t need to look for other systems- focusing on H6 chain is enough.

and it’s already good enough to have your correct FCI matrix, and I believe both we and Garnet want to step to the next stage ASAP. Just bear in mind what can be optimized - upper-diagonal, sparsity, etc. Garnet can ask you these and he care about these, which I do appreciate.