Spurious bader charges

mcorralv · Post by **mcorralv** » Tue Jan 17, 2006 2:03 pm

Hello,

I computed a CHGCAR for the ethylene molecule with VASP. I think the coordinates in the POSCAR file are fine as far as symmetry is concerned because VASP assigns the D2h point group to the molecule. However, when I compute bader charges with the CHGCAR in the output I get two different electron populations for both carbon atoms: 4.1566 and 3.8916.

The electron populations at the hydrogen atoms seem to be fine, only 1E-3 e- difference between them.

I don't understand why the bader charges in the carbon atoms are so different. The molecule is placed at the edge of the periodic cell, one carbon atom is one one side of the box while its neighbourg is at the other end of the box. Could that be the explanation? Does the program take into account the periodicity of the cell?

Thank you

MCV

Post by **graeme** » Tue Jan 17, 2006 3:22 pm

The code does take into account boundary conditions.

My first guess is that your fft grid needs to be finer. If the carbon atoms are lying in different positions with respect to the grid, the integrated charge for each atom will be different.

I suggest that you look in the OUTCAR or the CHGCAR to find the values of NGX, NGY, and NGZ. This is the size of your charge density grid. Then do a single ionic step run specifying NGXF, NGYF, and NGZF in your INCAR to be twice the original values. If this worked, you should get a CHGCAR file which is 8 times as large. Now rerun the bader analysis to see if the C charges are more accurate.

I would also appreciate it if you post back with your new charge numbers. This is a pretty new code, and we want to make sure it's working properly.

mcorralv · Post by **mcorralv** » Wed Jan 18, 2006 10:24 am

Hello,

Thank you for your reply. I tried again with a 360x360x360 mesh in a 18x18x18 periodic cell. I get an error message: "Segmentation fault". I tried in three different computers.

Should I try in a smaller cell or a in a computer with larger memory?

Thank you.

MCV

Post by **graeme** » Wed Jan 18, 2006 3:35 pm

360^3 is a pretty large data set, but I'm still surprised that it's causing a segmentation fault. Perhaps you are right about the memory. If you made the CHGCAR available I could try it on a large memory machine. You could also try one of your other suggestions. Reducing the cell size and the grid size would give you a higher density of points without increasing the total memory. You could also try a 1.5x increase in the grid size over the default grid -- perhaps 250^3.

andri · Post by **andri** » Wed Jan 18, 2006 7:23 pm

I have never gone over 200^3. Usually something of the order 150^3 is quite sufficient.

mcorralv · Post by **mcorralv** » Wed Jan 18, 2006 8:31 pm

My initial CHGCAR had a 180^3 grid.

Post by **graeme** » Wed Jan 18, 2006 8:46 pm

The number of necessary grid point is related to the box size. It is really the density of grid points that matter. An 18x18x18 Ang box is quite large, so a large number of grid points could be appropriate. However, if you are looking at an isolated molecule, this is a very large box. You could use a 10x10x10 Ang box with a 200x200x200 (or smaller) grid and still have a much higher mesh density than you have now.

Post by **graeme** » Fri Jan 20, 2006 5:44 am

Well, I think were stuck on this for now. Andri and I have tried a range of grid sizes and different INCAR parameters. I don't see a difference as large as you get. My biggest difference is 0.07e between C atoms, or 2%. Still, this is annoying. Which C atom has the larger charge also seems to depend on the grid size. The numbers don't seem to converge for grids as large as 224^3, so it doesn't look like an issue with grid density.

Anyways, we'll keep trying, but it was not the easy problem to fix that I was expecting. I don't think this has to do with the Bader partitioning -- you can see the same asymmetry in the Voronoi numbers. I think the asymmetry is really in the CHGCAR.

Do post back if you learn any more about this.

mcorralv · Post by **mcorralv** » Sat Jan 21, 2006 4:35 pm

I computed again the ethylene molecule with a different orientation of this
molecule in the cell: I turnt the molecule inside the cell so that none of its
symmetry planes were aligned with the cell's walls. Bader charges seem now to
be more reasonable.

# X Y Z VORONOI BADER % MIN DIST
---------------------------------------------------------------------------
1 7.2997 4.0650 5.3033 1.4750 0.9519 7.9329 0.3500
2 5.5619 4.9022 6.8637 1.4745 0.9808 8.1737 0.3239
3 5.9773 3.3754 4.2005 1.4745 0.9667 8.0561 0.3091
4 4.2395 4.2126 5.7609 1.4746 0.9615 8.0122 0.3609
5 6.2360 3.9141 5.1133 3.0505 4.0639 33.8661 0.5450
6 5.3032 4.3635 5.9509 3.0508 4.0751 33.9591 0.5551

This is the POSCAR file I used in the latter calculation:

header
1.00000000
10.000000 0.000000 0.000000
0.000000 10.000000 0.000000
0.000000 0.000000 10.000000
4 2
Cartesian
7.299694 4.065019 5.303280
5.561850 4.902170 6.863681
5.977339 3.375406 4.200527
4.239495 4.212558 5.760928
6.236021 3.914102 5.113301
5.303167 4.363475 5.950906

Kind regards,

Manuel

egillsk · Post by **egillsk** » Thu Feb 23, 2006 10:39 am

This is really strange. Is it possible that there is some bug in the boundery conditions of the algorithm?

Post by **graeme** » Thu Feb 23, 2006 4:46 pm

The molecule is in the center of the box, though, so isn't it unlikely that the boundaries are causing a problem?

egillsk · Post by **egillsk** » Fri Feb 24, 2006 1:04 pm

Well, initially when mcorralv calculated the molecule which was lying through the wall he got 4.1566 and 3.8916 electrons for the two different C atoms. When he calculated it when it was inside the cell, don't going through boundary conditions, then he got 4.0639 and 4.0751 el for each C atoms, which is off-course not as asymmtric. However, he probably had a finer grid in the later calculations which could be (and is probably) the reason. I just thought this was strange when I was reading through this.

andri · Post by **andri** » Fri Feb 24, 2006 7:03 pm

It just seems to me that if the boundary conditions were at fault here then there should also be a larger discrepancy between the hydrogen atoms initially when the molecule crosses the box boundary.

But there is something curious going on here that needs further investigation.

Post by **graeme** » Thu Sep 28, 2006 5:34 pm

This problem has finally been solved, thanks to Ed Sanville. He pointed out that our original algorithm biased the Bader surfaces to be along the grid directions. This bias leads to the deviations in charges for different orientations of the ethylene molecule. A modification to our algorithm, which is now the default in the bader code in cvs, gives consistent charges for the C atoms in ethylene (of ~4.05 e), independent of the orientation of the molecule.