voronoi construction problem

Dead lammps comunity,

I am learning how to perform voronoi construction in lammps and I am very confused. I add voronoi calculations to a simple LJ script. The density and temperature of LJ are deeply inside fcc phase, so almost a perfect crystal is expected. For the initial configuration (fcc crystal) I obtain that each atom has 12 nearest neighbors. However, already after the first step I see a distribution of nearest neighbors up to 17. Obviously, after just one step the system is stil a perfect crystal.

Does anyone know what is the reason of this problem? Any help is appreciated.

My input file is:

units lj
atom_style atomic

lattice fcc 1.2
region box block 0 6 0 6 0 6
create_box 1 box
create_atoms 1 box
mass 1 1.0

velocity all create 0.1 87287

pair_style lj/cut 2.5
pair_coeff 1 1 1.0 1.0 2.5

neighbor 0.3 bin
neigh_modify every 20 delay 0 check no

fix 3 all nvt temp 0.1 0.1 0.01

compute v1 all voronoi/atom
dump d1 all custom 1 dump.voro id type xu yu zu c_v1[1] c_v1[2]
compute r0 all reduce sum c_v1[1]
thermo_style custom c_r0
variable t1 equal c_r0

thermo 1
run 1

The results for the first configuraion:

0 8.96281
0 8.96281
0 8.96281
ITEM: ATOMS id type xu yu zu c_v1[1] c_v1[2]
1 1 0.000967711 0.000149307 8.00097e-05 0.83334 13
2 1 0.745183 0.747681 -0.00119005 0.833569 14
3 1 0.746222 -0.00245589 0.748209 0.833566 14
4 1 0.000175055 0.748187 0.745746 0.833895 15
5 1 1.49361 0.00101097 -7.97191e-06 0.833365 14
6 1 2.24308 0.747475 0.00244054 0.833676 14
7 1 2.23903 -0.00029294 0.744655 0.832089 13
8 1 1.49283 0.744114 0.744172 0.833446 15
9 1 2.98685 -0.00263345 -0.000921374 0.833923 14
10 1 3.73515 0.745935 0.00173668 0.835065 16
11 1 3.73268 0.00257329 0.748149 0.832109 12
12 1 2.98501 0.746174 0.747306 0.832129 12
13 1 4.47929 0.00159462 0.000950542 0.833319 13
14 1 5.22726 0.745646 0.00197512 0.833404 14
15 1 5.23024 -0.0018725 0.744675 0.834316 15
16 1 4.48155 0.746834 0.746483 0.832556 12
17 1 5.97261 -0.00257768 -0.00262245 0.832311 14
18 1 6.72357 0.747479 0.00157167 0.834482 13
19 1 6.72021 0.00220925 0.747646 0.833104 15

Well, if all of a sudden a particle has 17 nearest neighbours, I don’t think it is that obvious that the system should still be a perfect crystal. Have you visualised your atom positions? What happens to them? What does your potential energy look like?

Dear Stefan,

I am doing the simplest test possible. The system has done just one step with dt=0.005 LJ units. Nothing can happen in one step. Here is the log.file:

LAMMPS (9 Dec 2014)
Lattice spacing in x,y,z = 1.4938 1.4938 1.4938
Created orthogonal box = (0 0 0) to (8.96281 8.96281 8.96281)
1 by 1 by 1 MPI processor grid
Created 864 atoms
Setting up run …
Memory usage per processor = 3.61991 Mbytes
Step Temp E_pair E_mol TotEng Press
0 0.1 -7.6089166 0 -7.4590903 12.044613
1 0.09845862 -7.6064391 0 -7.4589221 12.057639
Loop time of 0.0635428 on 1 procs for 1 steps with 864 atoms

Pair time () = 0.00173092 (2.72402) Neigh time () = 0 (0)
Comm time () = 2.12193e-05 (0.0333936) Outpt time () = 0.0617428 (97.1672)
Other time (%) = 4.79221e-05 (0.075417)

Nlocal: 864 ave 864 max 864 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 2565 ave 2565 max 2565 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 57888 ave 57888 max 57888 min
Histogram: 1 0 0 0 0 0 0 0 0 0

Total # of neighbors = 57888
Ave neighs/atom = 67
Neighbor list builds = 0
Dangerous builds = 0

The final configuration is indistinguishable from the initial one. But the number of nearest neighbors definitely goes wrong.

The dump file is in attachment. It contains coordinates and both output values of compute_voro.

12.05.2015, 11:24, “Stefan Paquay” <[email protected]…24…>:

dump.voro (71.1 KB)

This is very strange. I figured you might have bad dynamics, but indeed, the atoms barely move but their coordination number is all over the place.
I am not at all an expert in crystals, but is the coordination in an fcc degenerate somehow? If I change the lattice to bcc, the coordination number of most particles is 14 and when I increase the temperature I see some having 13 and some 15, so it is as if it works correctly for a bcc lattice…

Just a follow-up, Ovito (a visualisation/analysis program) can perform common neighbor analysis, and with an adaptive variable cutoff it correctly identifies your system as FCC, so I think relying on Voronoi analysis to identify an FCC lattice is the real problem, rather than the Voronoi-compue itself.

I am not familiar with compute voronoi but is there a cutoff of some sort you can specify? 17, or 18, 19 does not seem odd to me for a fcc since the number of second nearest neighbors in fcc is 6. It could be possible that 1st and 2nd NN are mixed together…


You can specify a face threshold. I guess you accidentally count second nearest neighbours, but this will lead to a very small contact area, so if you would set this value to something slightly above 0 you might get rid of the problem.

Dear Stefan and Ray,

Thank you for your comments. However, they do not fix the problem.

First of all, the point is not to employ Voronoi construction to the fcc lattice. fcc is the simplest examle with known properties, so one can easily test the method with it. If we do not get correct results to fcc then I’m suspucious that the algorythm is correct.

Second, may be the second neighbors are involved in the calculations. However, then why do I get 12 nearest neighbors (the first shell only) for the initial configuration and up to 17 for the second one? Is the cut-off distance changed automatically after the first step? It looks very fishy. Moreover, by definition of the Voronoi construction only the first shell of the nearest neighbors should be counted.

The number of nearest neighbors in the bcc lattice is 8. The number 14 corresponds to the sum of the first and the second neighbors shells. So it also looks strange that we get 14 neighbors for the bcc structure.

If we try to compute the coordination of the simple cubic structure then we again meet the same problem as with fcc: in the initial configuration the coordination number of each particle is 6 (correct) while after the first tiny step which does not actually change the configuration the coordination number is up to 19.

Having that in 3 simplest test cases the command compute_voronoi gives definitely wrong results I come to one of two possibilities: (1) there is a bug in the procedure (2) something is wrong with my implementation.

Concerning the implementation - I just took it from the test case included into lammps distribution. Is there a bug in this test example? Or one should change something in order to empoy it to a different system?

Do you have any comments about this?

Kind regards,

12.05.2015, 17:05, “Stefan Paquay” <[email protected]…24…>:

Moreover, by definition of the Voronoi construction only the first shell of the nearest neighbors should be counted.

That’s not always the case. A Voronoi tesselation is a geometric entity.

For a (even slightly) disordered lattice, there is no strict definition of

“nearest neighbors”. Think of a 2d square lattice. There are 4 nearest

neighbors. But if every atom displaces randomly by epsilon, then the

Voronoi cells for each atom may have anywhere from 4-8 “neighbors”,

which include neighbors from what you are thinking of as 1st and 2nd

neighbor shells in a perfect lattice.

Unless you can show there is a bug in the way LAMMPS calculates

the Voronoi tesselation, this is not a LAMMPS issue. It’s an issue

with how you are interpreting the (correct) results.


Indeed, I think the issue is exactly what Steve is describing. With Voronoi tesselation, you just draw lines between atoms, and halfway that line you place a fictious, infinite plane with as normal the distance vector between the atoms. After this, each atom will be enclosed in some polyhedron, and the number of faces of that polyhedron is what the compute gives you. There is no cut-off for Voronoi tesselation, it is just the smallest enclosing polyhedron around each atom.

After tesselating, the number of faces for each polyhedron is counted and this is what the fix gives you, so if, due to some slight displacement there appears an additional face for that atom, it is counted as having an extra neighbour. I guess fcc and cubic lattices are susceptible to this, while a bcc is not, based on my observations. Just try yourself, replace “fcc” with “bcc”, and then it is correct. A bcc Voronoi cell has 14 faces, which consist of first and second nearest neighbours.

The issue is that at finite temperatures the number of faces of the voronoi cell in an fcc lattice fluctuates heavily because these small displacements create extra faces. One way to work around this issue of “sharp” Voronoi cells is by placing a threshold on how large a face needs to be before it is counted as an additional neighbour (this is what the keyword face_threshold sets).

Not sure what you are after with the Voronoi tool but if by some reason your goal is to find a descriptor to locally classify a lattice structure and/or deviations from it take a look at:


Its an angular-based descriptor less susceptible to temperature effects. A bit limited on the classification range but maybe that’s within your similarity error bars.