Issue with Periodicity for Bonds

Hi Axel,

I’m encountering a problem with the periodicity of a 2D perovskite, and I’ve attached the relevant files in cif, cfg, and pdb formats to this message. Specifically, I’m trying to use topotools to extract topological information such as bonds, angles, and dihedrals, but I’ve noticed that some bonds are missing at periodic boundaries along the Y axis. I posted about this issue on the VMD list, but haven’t received any helpful responses yet. So, I was wondering if you might have any suggestions for addressing this problem within LAMMPS itself.

Thanks a lot for taking the time to consider this matter.

Best regards,

Reza.

n2.zip (77.3 KB)

Hi,

I am not Axel, but I can give you a solution that is not LAMMPS nor VMD related. You can very easily build a python (or whatever language you prefer) code yourself to create the connectivity of your system. This is an easy thing to do, especially if your atoms are already in optimal position relative to one another, which seems to be your case. Also you will have full control on which bonds/angles/dihedrals are being input (and also if it is being done correctly) since you are the person writting the code.

You can set up a distance criteria to create the bonds - for example if two atoms of a given type are sitting at a distance smaller than a threshold (dictated by the bond length), you can create an array containing their atom-IDs to later declare it in the connectivity section of the MD software you are using. Of course you need to take into account the fact that meant-to-be-bonded atoms may be sitting at different edges of the simulation box as your system is periodic. I can give you a small, very step-by-step example on how you can check for that within python and then you can optimize it and build the rest of the code.

# atom_1 and atom_2 are 1x2 vectors containing position in x and y direction.
# i and j are also 1x2 "generators" of your simulation box. In your case i[0,0] carries the lx value and j[0,1] caries the ly value.
dist_x = atom_1[0,0]-atom_2[0,0]
dist_y = atom_1[0,1]-atom_2[0,1]
if (abs(dist_x) > abs(((mod_i/2)))) & (dist_x > 0) :
    dist_x = dist_x - i[0,0]
if (abs(dist_x) > abs(((mod_i/2)))) & (dist_x < 0):
    dist_x = dist_x + i[0,0]
if (abs(dist_y) > abs(((mod_j/2)))) & (dist_y > 0) :
    dist_y = dist_y - j[0,1]
if (abs(dist_y) > abs(((mod_j/2)))) & (dist_y < 0):
    dist_y = dist_y + j[0,1]
dist_total=((dist_x)**2 + (dist_y)**2 + (dist_z)**2)**(1/2)

You should be able to tell how it works by drawing in a piece of paper a scenario in which atom_i and atom_j are sitting at edges of the simulation box. And then afterwards you can easily build angles and etc once you have the bond section. Maybe not as fancy as using other software to build the connectivity but it works quite well!

1 Like