Top Page > Browsing
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.
メンテ
Page: [1]

Re: Wrong calculation of chemical potential in some cases ( No.1 )
Date: 2020/12/10 18:38
Name: Naoya Yamaguchi

Hi,

I think that your point is right.
You can find the following in "Band_DFT_Col.c" in OpenMX 3.9.2.

lumos = (double)n*0.200;
if (lumos<60.0) lumos = 400.0;
MaxN = (TZ-system_charge)/2 + (int)lumos;
if (n<MaxN) MaxN = n;

I guess that the developers prepared "Band_DFT_NonCol.c" from "Band_DFT_Col.c", and they forgot to modify the line of the first determination of "MaxN".
And, the line is clearly wrong, so if I checked that in OpenMX 3.8.5, the corresponding part is

lumos = (double)2*n*0.100;
if (lumos<60.0) lumos = 400.0;
MaxN = (TZ-system_charge) + (int)lumos;
if ( (2*n)<MaxN ) MaxN = 2*n;

The most important point is that "/2" in the line must be deleted, and the first determination of "lumos" should be modified, if necessary.

Owing to your post, I noticed that this issue was the same as Nos. 8-12 in the following thread.
http://www.openmx-square.org/forum/patio.cgi?mode=view&no=2671

Regards,
Naoya Yamaguchi
メンテ

Page: [1]

Thread Title (must) Move the thread to the top
Your Name (must)
E-Mail (must)
URL
Password (used in modification of the submitted text)
Comment (must)

   Save Cookie