lmhale99 / atomman

Atomistic Manipulation Toolkit
Other
34 stars 40 forks source link

Neighborlist returns `10` instead of `12` as coordination number for hcp material. #11

Closed mastricker closed 4 years ago

mastricker commented 4 years ago

Hello,

during the analysis of a dislocation system with a hcp crystal structure I encountered an issue with the neighorlist, which is needed e.g. for the p_vectors for the Nye tensor evaluation.

hcp_system in my context is an orthogonal box with the basal plane in the x-y plane. The specific lattice parameter in my case are: a = 3.1790434213597902, c = 5.100201896236235

Following code should return 12 as the number of p_vectors because it is the number of nearest neighbors in a hcp system but it returns 10:

The cutoff reflects the distance of atoms in the basal planes - diagonally across to the other planes is shorter than a for the nearest neighbors.

cutoff = 1.01 * a

# load reference system
hcp_system = am.load_atom_data("atom_ref.dat")
hcp_system.pbc = (True, True, True)

neighbors0 = hcp_system.neighborlist(cutoff=cutoff)

all_p_vectors = am.defect.nye_tensor_p(hcp_system, neighbors=neighbors0)
p_vectors = all_p_vectors[0]

unique_p_lens = set([len(p) for p in all_p_vectors])

print("Unique p vector lengths (should be 12 for hcp?) = {}".format(unique_p_lens))

When I analyse the same system with the ASE (atomic simulation environment) neighborlist algorithm, I get the correct 12 number of nearest neighbors for the given cutoff.

In short: if this is not correct, I assume all following analysis is not correct, too, because 2 p_vectors are missing.

Thank you. Let me know if you need any more info.

Best, Markus

PS: I'm using atomman-1.3.1 with python3.7.

lmhale99 commented 4 years ago

Hi Markus,

I just ran the test below and get 12 neighbors:

a = 3.1790434213597902
c = 5.100201896236235

# Build primitive HCP ucell
atoms = am.Atoms(pos = [[0.0, 0.0, 0.0], [1/3, 2/3, 1/2]])
box = am.Box.hexagonal(a=a, c=c)
pcell = am.System(atoms=atoms, box=box, scale=True)

# Convert to conventional HCP ucell
ucell = pcell.rotate([[1, 0, 0],[1,2, 0],[0, 0, 1]])
print(ucell)

# Increase size and shift
system = ucell.supersize((3), (3), (3))
system.atoms.pos -= np.array([0.0, 0.0, 1.0])
system.wrap()

# Identify neighbors and cutoff
neighbors = system.neighborlist(cutoff=1.01*a)
print(neighbors.coord.mean())

The two possibilities that I can think of for your results are that

If I use smaller values in supersize, I can get coordinations between 3-12, but larger values always give coordinations of 12. Try replicating your system along any short dimensions.

In my tests, I did see an actual bug that neighbors near the x and z boundaries were occasionally missed. This meant that the mean coordination was coming out to be around 11.9 rather than 12, so it shouldn't be related to what you are seeing. The shift that I apply after the supersize seems to avoid the bug and give the correct mean of 12.

Let me know if this solves the problem, Lucas

mastricker commented 4 years ago

Hi Lucas,

thanks for the quick response. I ran it with the suggested replication in the shortest direction and indeed, the neighborlist is now showing the correct 12 neighbors.

Thanks again! Best, Markus