The problem is my box is not orthorhombic (base angles different than 90 degrees). So, I am getting wrong results with the formula,

def distance(a, b):
    dx = abs(a[0] - b[0])
    x = min(dx, abs(A - dx))
    dy = abs(a[1] - b[1])
    y = min(dy, abs(B - dy))
    dz = abs(a[2] - b[2])
    z = min(dz, abs(C - dz))
    return sp.sqrt(x**2 + y**2 + z**2)

Any suggestion on how to get the proper distance between any atom and its neighbours, including the periodic images of the neighbours. Basically, it makes the results of RDF and CN unphysical.

enter image description here

  • 1
    $\begingroup$ Have you tried using the Ovito tool? There are built in functions to evaluate the RDF and CN from the menu. Alternatively, it might be possible to use a tool such as Vesta and convert the unit cell to a rectangular one. $\endgroup$
    – PBH
    Dec 1 '21 at 1:33
  • $\begingroup$ @PBH why not write an answer? Comments are temporary and can be removed permanently without any notification. $\endgroup$ Dec 1 '21 at 16:48

Sticking with ase you can use the built-in get_distance function. It supports non-orthogonal unit cells.

An example to get the distance between atom number a and atom number b is:

import ase.io
atom = ase.io.read('your_file.xyz')
distance_a_to_b = atom.get_distance(a, b)

Additionally if you want to compute the Radial Distribution Function (RDF) and Coordination Number (CN), asap3 is a python package building on top of ase that natively does it.

  • $\begingroup$ I used the asap3 for the RDF calculation but the results from it doesn't give correct CN. Also, surprisingly, the RDF for the same file is not the same for OVITO and ASAP3, with same bin size, cutoff etc. $\endgroup$
    – Saha_1994
    Dec 2 '21 at 3:49
  • $\begingroup$ Have you checked that the cell was properly read by ase (atom.get_cell()) and periodic boundary conditions set up (atom.get_pbc()) ? $\endgroup$
    – Hebo
    Dec 2 '21 at 14:54
  • $\begingroup$ Hi, @hebo, upon checking I found that the cell is read properly bu ASE along with the pbc. Even more baffling. $\endgroup$
    – Saha_1994
    Dec 14 '21 at 7:14

Have you tried using the Ovito where there are built in functions to evaluate both RDF and CN of a system?

Additionally, it might be possible to use Vesta to convert the unit cell to a rectangular one (see this video). That way you can use your code to do the analysis.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.