Wrong calculation of chemical potential in some cases |
- Date: 2020/12/03 13:56
- Name: Seungjin Kang
<gsj37@snu.ac.kr>
- I have encountered some cases where the chemical potential is not calculated correctly.
- 1. Input file
I used the the input file for Fe-bcc in openmx work directory, openmx/work/large_example/FeBCC.dat
Then I changed some parameters as follows.
scf.SpinPolarization nc
scf.SpinOrbit.Coupling on
scf.maxIter 1
This turns on non-collinear calculation, and scf.maxIter is just set to be 1 since the only purpose is to check whether chemical potential is correctly calculated or not.
- 2. Test case
This system contains 16 Fe atoms.
For test calculation, I varied the number of Fe atoms by deleting some lines in Atoms.SpeciesAndCoordinates section.
With 15 Fe system, for example, I just deleted
[16 Fe 0.75 0.75 0.75 8.0 6.0]
line in Atoms.SpeciesAndCoordinates, and for 14 Fe system, deleted
[15 Fe 0.75 0.75 0.25 8.0 6.0]
[16 Fe 0.75 0.75 0.75 8.0 6.0]
and so on. Here's the result of the chemical potential value from the output file.
[11 Fe atoms case]
Chemical Potential (Hatree) = -0.15186311662247
total= 154.00000 ideal(neutral)= 154.00000
[12 Fe atoms case]
Chemical Potential (Hatree) = 20.00000000000000
total= 146.00000 ideal(neutral)= 168.00000
[13 Fe atoms case]
Chemical Potential (Hatree) = 20.00000000000000
total= 158.00000 ideal(neutral)= 182.00000
[16 Fe atoms case]
Chemical Potential (Hatree) = 20.00000000000000
total= 195.00000 ideal(neutral)= 224.00000
So, up to 11 Fe atoms, the chemical potential is calculated correctly, so the number of electrons is also correctly calculated.
But if we add more Fe atoms, we see that the chemical potential and the number of electrons is not calculated correctly.
In fact, the number 20.0000 for chemical potential is the default ChemP_MAX value.
- 3. Origin of the Problem
I guess that the origin of this problem is the wrong estimation of the number of state to be solved in the calculation.
For non-collinear case, this is determined between the line 248~251 in Band_DFT_Noncol.c file.
[Band_DFT_NonCol.c]
248 lumos = (double)n2*0.200;
249 if (lumos<60.0) lumos = 400.0;
250 MaxN = (TZ-system_charge)/2 + (int)lumos;
251 if (n2<MaxN) MaxN = n2;
So, MaxN defines the number of states to be solved.
And I found that For the Fe-bcc case, since our basis set is s2p2d1 = 13*2=26, the n2 (Total number of basis set) changes as follows.
Number of Fe atoms - n2 - lumos
11 - 286 - 57.2
12 - 312 - 62.4
13 - 338 - 67.6
After 12 atoms, the value of lumos exceeds 60 and line 249 is ignored. Then MaxN is set to be
Number of Fe atoms - MaxN - Total number of Electrons
11 - 286 - 154
12 - 146 - 168
13 - 158 - 182
Since MaxN is smaller than the total number of electrons after 12 atoms, it is apparent that the chemical potential cannot be calculated correctly.
The easiest way to resolve this issue is to set MaxN = n2, as
MaxN = n2;
after line 251.
I am not sure what would be the optimal way to choose MaxN, maybe we can just increase lumos value.
But I think this issue could also occur for other systems (supercell calculation, etc.), so proper handling of MaxN is required.
 | |
This thread is locked.Only browsing is available.
Re: Wrong calculation of chemical potential in some cases