Re: overlap between k1 and k2 ( No.1 ) |
- Date: 2022/05/10 10:42
- Name: Naoya Yamaguchi
- Dear Wei,
The `Sop` is eq. 11 of http://www.openmx-square.org/tech_notes/tech11-1_1.pdf
The braket of the periodic parts of the Bloch functions at two different k-points should not be 0.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.2 ) |
- Date: 2022/05/10 18:13
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
Thanks so much for your helps.
Really helps me a lot.
Wei
|
Re: overlap between k1 and k2 ( No.3 ) |
- Date: 2022/05/11 05:02
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
I just added few lines after call Overlap_k1k2 to make it print the Sop.
for example,
for (spin=0; spin<spinsize; spin++){ for (i=0; i<fsize3; i++){ for (j=0; j<fsize3; j++){ fprintf(f1, " %.2e %.2e", Sop[spin][i][j].r, Sop[spin][i][j].i); } fprintf(f1, " \n "); } } //spinsize
Just found most of the elements are zero, few are non-zero. Is this because the Sop is in AO basis, and overlap integral of some AOs are very small?
1st line: all zero
2nd line: -2.21e-04 7.88e-05 -1.29e-08 4.16e-09 1.25e-11 9.09e-11 9.03e-04 -6.23e-03 -4.56e-11 -1.64e-11 -1.20e-01 8.15e-02 4.55e-07 -1.66e-07 1.14e-01 -8.47e-03 5.30e-10 7.43e-11 1.94e-11 1.78e-11 -1.23e-03 4.57e-04 -1.87e-09 7.14e-10 2.20e-03 -3.47e-04 -2.36e-12 -5.14e-11 1.90e-03 7.53e-05 -4.11e-11 2.60e-11 2.01e-03 -3.56e-04 4.75e-11 -2.95e-11 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00 0.00e+00
3rd line: similar to 2nd line
4th line: similar to 2nd line
...
thank you very much.
Wei
|
Re: overlap between k1 and k2 ( No.4 ) |
- Date: 2022/05/11 12:28
- Name: Naoya Yamaguchi
- Dear Wei,
The returning `Sop` depends on the arguments of `Overlap_k1k2`, and I can't say anything unless you provide them.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.5 ) |
- Date: 2022/05/11 18:25
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
My system is CsPbI3 unit cell with SOC included. I compiled polB.c, and providing the following inputs:
Specify the number of grids to discretize reciprocal a-, b-, and c-vectors (e.g 2 4 3) 2 4 3
Specify the direction of polarization as reciprocal a-, b-, and c-vectors (e.g 0 0 1 ) 0 0 1
Just for a quick test, and add the lines in No.3 for printing Sop. The overlap I gave in No.3 is between two different kpoints at k1 0.500000 0.750000 0.000000 and k2 0.500000 0.750000 0.3333.
Thank you. Wei
|
Re: overlap between k1 and k2 ( No.6 ) |
- Date: 2022/05/12 02:29
- Name: Naoya Yamaguchi
- Dear Wei,
I wanted to know what was an argument of `diag_flag` for the function `Overlap_k1k2` because the `diag_flag` should be set as appropriate in calling it.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.7 ) |
- Date: 2022/05/12 05:38
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
Q1: The dia_flag was set internally. I did no changes.
There are three loops for printing Sop: diag_flag: 0 overlap between k1 0.500000 0.750000 0.000000 k2 0.500000 0.750000 0.333333 do print Sop[spin][i][j] diag_flag: 2 overlap between k1 0.500000 0.750000 0.333333 k2 0.500000 0.750000 0.666667 do print Sop[spin][i][j] diag_flag: 3 overlap between k1 0.500000 0.750000 0.666667 k2 0.500000 0.750000 1.000000 do print Sop[spin][i][j]
each output is similar to case in No. 3, even for diagonal_flag=0. Namely, most are zeros, few are non zero.
Q2: btw, another question is associated with the H and S with non-collinear spin-polarized including SOC. As I see from Hamiltonian_Band_NC function, H is with the following form: for (i=1; i<=NUM; i++){ for (j=1; j<=NUM; j++){ H[i ][j ].r = H11r[i][j]; H[i ][j ].i = H11i[i][j]; H[i+NUM][j+NUM].r = H22r[i][j]; H[i+NUM][j+NUM].i = H22i[i][j]; H[i ][j+NUM].r = H12r[i][j]; H[i ][j+NUM].i = H12i[i][j]; H[j+NUM][i ].r = H[i][j+NUM].r; H[j+NUM][i ].i = -H[i][j+NUM].i; } }
here NUM is the number of AO basis. So the H matrix will be 2NUMx2NUM in size. Also, the H matrix can be regarded as a block matrix in the following form: Haa Hab Hba Hbb
S is still a NUMxNUM matrix.
We want to use subroutines of python to diagonalize H and obtain the eigenvalues and eigenfunction in MO basis through H * C = S * C * E However, now H and S have different size, and cannot do matrix multiplication. Do you have a solution here?
Generally, H with SOC can adopt the spinor form, not sure if it can be applied here. We expect the eigenfunction C is also in the spinor form.
Thank you very much.
Wei
|
Re: overlap between k1 and k2 ( No.8 ) |
- Date: 2022/05/12 12:25
- Name: Naoya Yamaguchi
- Dear Wei,
A1: In No. 3, you meant that Sop[0][i=0][j=*] is zero, didn't you? It is natural because the meaningful indices of i and j are from 1 in the released version in Ver. 3.9. You can find it from the bottom of the `Overlap_k1k2`.
A2: S 0 0 S
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.9 ) |
- Date: 2022/05/12 19:41
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
thank you.
In No.3 I meant: Sop[0][i=0][j=*] all zero Sop[0][I=1][j=0:k] non-zero Sop[0][I=1][j=k+1:fsize3] zero Sop[0][I=2][j=0:k] non-zero Sop[0][I=2][j=k+1:fsize3] zero ... Sop[0][I=n][j=*] all zero again
This happens in both non-collinear and collinear calculations.
thank you very much.
Wei
|
Re: overlap between k1 and k2 ( No.10 ) |
- Date: 2022/05/12 18:39
- Name: Naoya Yamaguchi
- Dear Wei,
The `Sop[0][1][0]` should also be zero, but your output was not true for it. In `polB.c`, `Sop[0][0][*]` or `Sop[0][*][0]` are not used in the calculations, and I guess that your additional modification makes the non-zero values. I can't say more without the entire code you modified.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.11 ) |
- Date: 2022/05/12 19:44
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
Thank you. I will double-check.
Btw, just to confirm with you: I have obtained the identical SOC eigenvalues as openmx outputs by diagonalizing H . Then the eigenfunction (or wave function) is also in the following form, right? WF_aa WF_ab WF_ba WF_bb
Since H or S is in: Haa Hab Hba Hbb
S 0 0 S
Thank you.
|
Re: overlap between k1 and k2 ( No.12 ) |
- Date: 2022/05/12 21:13
- Name: Naoya Yamaguchi
- Dear Wei,
If you mean C, it is not correct. The C stores the LCPAO expansion coefficients for each eigenvalue in each column. The two spin components appear in each column: the first half is alpha, and the reminder is beta.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.13 ) |
- Date: 2022/05/13 02:55
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
Thank you. The size of C matrix is 2N * 2N, here N is the number of AO basis. We want to obtain the spinor wave function by diagonalizing H.
Isn't the C matrix the spinor wave function? If it is, then I am a bit confused on the matrix elements. You mentioned each columns appear in the order of alp bet alp bet alp bet ... with each two columns (alp and bet) correspond to the spinor wave functions. here C is a square matrix, which means that the size of spinor wave function of corresponding MO is also 2*N, right? Please let me know if I am wrong. My understanding here is the wave function is in MO basis, instead of AO basis now.
If C is not the spinor wave function, could you please point me how to extract it. I understand polB or calB can probably get the spinor wave function at corresponding k (e.g., Wk1 or Wk2), but I afraid it maybe difficult for me to do further development based on polB or calB.
Best regards
Wei
|
Re: overlap between k1 and k2 ( No.14 ) |
- Date: 2022/05/13 21:05
- Name: Naoya Yamaguchi
- Dear Wei,
>Isn't the C matrix the spinor wave function? >If C is not the spinor wave function, could you please point me how to extract it.
The C stores c in eq. 1 of http://www.openmx-square.org/tech_notes/tech1-1_2.pdf
>You mentioned each columns appear in the order of alp bet alp bet alp bet
No, the first N columns are alpha, and the reminder are beta. It is clear from the relation of the product of matrices of HC.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.15 ) |
- Date: 2022/05/15 17:35
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
Thank you. Seems I get a bit confused previously. Again, thanks so much for pointing me the right way.
Actually we want the overlap between orbital i and j at different k, Psi_ij^k1k2, at different time t and t+dt, considering we have different snapshots along MD.
Q1: For the simple case, dt=0, and for the same snapshot: Psi_i^k1 = e^ik1 R C_i^k1 phi_i; Psi_j^k2 = e^ik2 R C_j^k2 phi_j;
then according to eq. 11 of tech11-1_1 note, the overlap:
Psi_ij^k1k2 = <Psi_i^k1 | Psi_j^k2> = C^k1^* C^k2 e^-ik1R e^-i(k2-k1)R {<phi_i|phi_j> - i(k2-k1)<phi_i|r|phi_j>}
Here, C matrix is eigenvector, LCAO matrix, and can be obtained from diagonalizing H, the exponent term is known. Here <phi_i|phi_j> is overlap S printed in HS.out file, right?
Q2: Now extend to the time overlap:
Psi_ij^k1k2 = <Psi_i^k1(t) | Psi_j^k2(t+dt)> = C(t)^k1^* C(t+dt)^k2 e^-ik1R e^-i(k2-k1)R {<phi_i(t)|phi_j(t+dt)> - i(k2-k1)<phi_i(t)|r|phi_j(t+dt)>}
Here, C(t) and C(t+dt) terms can be obtained from diagonalizing H at t and t+dt. exponents term is also know. But here I am not sure <phi_i(t)|phi_j(t+dt)>, evaluation of this term requires we read the AO wave function at different time, and do the AO overlap. Could you please give me some hints on this?
In some ab initio codes which uses the Gaussian type basis set, it generates the molden file, which includes the eigenvector, exponent and contractions coefficients of GTO. The exponent and contractions coefficients at different time can be used to compute the AO overlap at t and t+dt, but not sure about how Openmx code can print the similar infos for AO overlap calculation.
Best,
Wei
|
Re: overlap between k1 and k2 ( No.16 ) |
- Date: 2022/05/15 21:12
- Name: Naoya Yamaguchi
- Dear Wei,
A1: >Here <phi_i|phi_j> is overlap S printed in HS.out file, right?
Yes, it is stored in an scfout file.
However,
>Psi_ij^k1k2 = <Psi_i^k1 | Psi_j^k2>
"Psi" has to be the periodic part of a Bloch function.
A2: Generally, the overlap in the LHS of eq. 11 is gauge-dependent, and I think that adjusting the gauge between the two snapshots is difficult, because the snapshots are independent of one another, that is, the time evolution of the wave functions are not considered.
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.17 ) |
- Date: 2022/05/15 22:02
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
thanks for your reply.
Now clear for Q1.
For Q2, polB.c outputs WF at k1 and k2, the WF is the AO WF, right?
Just thinking if <phi_i(t)|phi_j(t+dt)> can be expressed as WF overlap. If yes, then probably we can store WF at t and t+dt, and then do the overlap. But not sure how to get the WF, though it has been done in polB.c. Now I can parse all variables stored in read_scfout.c using python.
if not, then probably need to think another way.
Our methodology development needs these time overlaps.
Best regards,
Wei
|
Re: overlap between k1 and k2 ( No.18 ) |
- Date: 2022/05/16 01:33
- Name: Naoya Yamaguchi
- Dear Wei,
While the "Psi" between diffrent k-points is difficult, "phi" is possible to do so, but the "overlap" (correlation) of "phi" at different times is not stored in scfout files. You can prepare it by modifying the source code of OpenMX by yourself instead of scfout files.
>For Q2, polB.c outputs WF at k1 and k2, the WF is the AO WF, right?
The `polB` evaluates the gauge-independent Berry phase (eq. 3), not the gauge-dependent quantity of eq. 11 itself: See also ref. 1 of http://www.openmx-square.org/tech_notes/tech11-1_1.pdf : "R. D. King-Smith, and D. Vanderbilt, Phys. Rev. B 47, 1651 (1993)".
Regards, Naoya Yamaguchi
|
Re: overlap between k1 and k2 ( No.19 ) |
- Date: 2022/05/16 05:11
- Name: Wei Li <liwei0099@gmail.com>
- Dear Naoya,
thanks for your quick feedback. For overlap at different time, sounds my thoughts on WF are wrong, I should move to "phi".
Q1: In order to store "phi", could you please point out which part of the code I need to look at?
I would like to do the following: 1) scf calculation at a snapshot (t), get "phi(t)", 2) scf calculation at another snapshot, get "phi(t+dt)", 3) do the overlap <phi_i(t)|phi_j(t+dt)>, and <phi_i(t)|r|phi_j(t+dt)>.
These are needed to evaluate <Psi_i^k1 (t) | Psi_j^k2 (t+dt)>, as in No. 15.
Q2: Evaluation of <phi_i(t)|r|phi_j(t+dt)> maybe tricky, since it needs the OLPpo term between different tilmestep ... may need to dig into the code how OLPpo was calculated using "phi". Could you please point me which part of the code calculates OLPpo and OLP?
Actually, we just need the size of phi(t) and phi(t+dt) can match to do matrix multiplication.
Q3: eq. 11 of tech11-1_1 tell us: <psi^k|psi^k'> = C^(k)* C^(k') e^ikRn e^-idk (-Rn) x {<phi_i|phi_j> - idk <phi_i|r|phi_j>}
As I check from polB.c script, the code first call setexpOLP() function, double phase=-(dkx*Gxyz[ct_AN][1]+dky*Gxyz[ct_AN][2]+dkz*Gxyz[ct_AN][3]); double co=cos(phase); double si=sin(phase); expOLP[ct_AN][h_AN][i1][j1].r=co*tmp3r-si*tmp3i; expOLP[ct_AN][h_AN][i1][j1].i=co*tmp3i+si*tmp3r;
So here means {<phi_i|phi_j> - idk <phi_i|r|phi_j>} multiplies exp^-i(dk)(-Rn)?
And then the code call expOLP_Band() function: kRn = k1*(double)l1 + k2*(double)l2 + k3*(double)l3; kRn-=dka*l1+dkb*l2+dkc*l3; T1[Anum+i][Bnum+j]+=t.r*co-t.i*si; T2[Anum+i][Bnum+j]+=t.i*co+t.r*si; T1[Bnum+j][Anum+i]+=t.r*co2+t.i*si2; T2[Bnum+j][Anum+i]+=t.i*co2-t.r*si2;
Does this mean expOLP multiplies exp^-i(dk)(-Rn) again?
Just thinking may need to multiply exp^-i(dk)(-Rn) first, then multiply exp^-ikRn. Perhaps I made something wrong... Just a bit more question, thanks so much for your helps.
Best regards,
Wei
|
Re: overlap between k1 and k2 ( No.20 ) |
- Date: 2022/05/16 15:48
- Name: Naoya Yamaguchi
- Dear Wei,
A1: To prepare the overlap between the "phi", that is, PAOs, the scf calculation is not necessary.
A2: Since the PAOs are expanded around atoms, you may take the overlap between those at two different times. >Could you please point me which part of the code calculates OLPpo and OLP?
You may modify the `SCF2File.c` if my memory is right, and you can dig into the source code of the OpenMX with the `grep` command.
A3: >So here means {<phi_i|phi_j> - idk <phi_i|r|phi_j>} multiplies exp^-i(dk)(-Rn)?
Yes, you should see eqs. 9 and 10 in the technical note.
The `expOLP` and `expOLP_Band` are based on a new implementation (the evaluation is not changed), and it follows Eq. 4 of https://arxiv.org/ftp/arxiv/papers/2203/2203.10441.pdf and T (below Eq. 4).
>T1[Anum+i][Bnum+j]+=t.r*co-t.i*si; >T2[Anum+i][Bnum+j]+=t.i*co+t.r*si; >T1[Bnum+j][Anum+i]+=t.r*co2+t.i*si2; >T2[Bnum+j][Anum+i]+=t.i*co2-t.r*si2;
These lines just gives a symmetrized T in terms of the k-points, and you can refer to the `sign==M` and `sign==P` cases as more simple examples.
Regards, Naoya Yamaguchi
|
|