Compute pair/local distance yields opposite sign than expected

LAMMPS LAMMPS Development
1

Hi,

I’m using the compute pair/local command to get the distances from a particular set of atoms of the system using the python interface. I’m using the 23Jun2022 (Update 1) release.

I set up two compute commands as:

lmp.command(“compute ids all property/local patom1 patom2”)
lmp.command(“compute dist all pair/local dx dy dz”)
lmp.command(“run 0”)

to get the id pairs of all interactions, and the distances associated to those interactions. lmp is the LAMMPS object.

Then I extract the information from those computes:

ids = gather_compute(“ids”,lmp,MPI.COMM_WORLD,2,np.int32)
dist = gather_compute(“dist”,lmp,MPI.COMM_WORLD,3)

gather_compute is a small function that I made to collect the information with MPI:

def gather_compute(computeid,lmp,comm,vec_size,datatype=np.float64):
    array = lmp.numpy.extract_compute(computeid,LMP_STYLE_LOCAL,LMP_TYPE_ARRAY).astype(datatype,copy=False)
    farray = array.flatten()
    array_sizes = np.array(comm.gather(farray.shape[0],root=0))
    if(master):
        recvbf = np.empty(array_sizes.sum(),dtype=datatype)
    else:
        recvbf = None
    comm.Gatherv(farray,(recvbf,array_sizes),root=0)
    if(master):
        return recvbf.reshape((int(array_sizes.sum()/vec_size),vec_size))
    else:
        return None

According to the LAMMPS documentation for both the latest stable release and the version I’m using the vectors in ‘dist’ should always be defined from the smallest id to the biggest id (“This value [dx,dy,dz] is always the distance from the atom of lower to the one with the higher id.”) but I find that the opposite is true. Outputting a couple i,j pairs of the previously defined ‘ids’ and their corresponding ‘dist’ vectors I get:

1539 562 [-3.94407078 8.67170072 0.98916273]
1542 562 [-2.12934197 8.51996973 2.73700122]

The positions of the atoms with these ids are:

562 -0.1641226737947051 8.984665968777039 -0.9482579825615813
1539 3.7799481076267973 0.31296524725565583 -1.9374207112908028
1542 1.9652192984527783 0.4646962400097415 -3.6852592019084223

As you can see for the 1539 562 pair if we did the vector from 562 (smaller id) to 1539 (bigger) we would get pos(1539)-pos(562) = 3.94407078, -8.67170072, -0.98916273, the opposite sing. This is easy to fix from my python script but I would like to know if this is the intended output.

Unless I’m misunderstanding the LAMMPS documentation it appears this is a bug? I have not found in the release notes if this has been fixed, but sorry if this is known/fixed in some later release.

Thanks.

2

It is difficult to assess the situation from the bits and pieces of your description.
Instead, it would be much easier, if you could construct a minimal set of test files that can be run and demonstrate the discrepancy directly. That would massively reduce any efforts to track down what the specific explanation is and whether this is incorrect behavior of the LAMMPS code or incorrect documentation or an incorrect interpretation of the data.

3

Hello Axel, apologies from my lack of clarity.

I have prepared a minimal example based on the simple.py example of the LAMMPS source code:
in.simple (376 Bytes)
simple.py (2.3 KB)
I have runned this example (with MPI) with mpirun -np 4 python simple.py in.simple
The code produces the output of a specific interaction pair that reproduces de problem, but I have checked and other pairs present the same problem (a PBC correction may be necessary when doing the difference in positions, I have chosen this pair so this is not necessary).

The output of the code on my machine is:

Interaction index i j distance
Interaction index:6656, i,j: 107 158, distance reported by compute pair/local: [-1.70686088 0.82582049 -0.91695824]
i (107) postion [4.26766699 6.67956782 5.800046 ]
j (158) postion [5.97452786 5.85374733 6.71700424]
pos(j)-pos(i) (vector from i:107 to j:158): [ 1.70686088 -0.82582049 0.91695824]

As you can see pair/local has a changed sign. Let me know if you need any more info, I hope this is clear enough.

Thank you.

4

Thanks for the files.

Yes, I can confirm that the output is contradicting the documentation.
The question is now whether to change the code or the documentation.

@sjplimp @athomps @stamoor what do you think?

5

I think the code is fine as is, and therefore should not be changed. I have not seen anything close to a convention in MD codes for preferring x1-x2 over x2-x1. Better to change the doc page to “This value is always the distance from the higher id atom to the lower id atom.”

1 Like