Fix bug in atom sorting with triclinic boxes

[stanmoore1]

Summary

Atom sorting for triclinic boxes is currently broken because atom positions are in 0...1 "lamda" coords but the bounding box and binsize are in regular box coords. This makes sorting a no-op for many cases since all atoms end up in the same bin. This PR fixes the issue by temporarily converting the atom positions back to box coords for sorting.

Alternatively the bounding box and binsize could be converted to lambda coords, but that seems more complex and would require many more changes to the code. We could consider this path though if other developers think it is better.

Related Issue(s)

#3740

Author(s)

Stan Moore (SNL)

Licensing

By submitting this pull request, I agree, that my contribution will be included in LAMMPS and redistributed under either the GNU General Public License version 2 (GPL v2) or the GNU Lesser General Public License version 2.1 (LGPL v2.1).

Backward Compatibility

Yes

[stanmoore1]

Looking at this more, it seems like the "correct" solution would be to bin in lamda coords (similar to PPPM) instead of using a bounding box which has extra volume for a skewed box. However, this would require changing the neigh list binning too so everything matches to maximize cache hits/data locality.

[stanmoore1]

Though lamda coords with a highly skewed grid could bin up far away atoms in the corners of the skewed bins. I think we should just merge this PR for now.

[athomps]

I did some testing on this, and verified that you can bin in lambda coords, as you suggest. The key is to define bbox in lamda coords, which are provided by domain->sublo_lamda, etc.:

  if (domain->triclinic) {
    //    domain->bbox(domain->sublo_lamda,domain->subhi_lamda,bboxlo,bboxhi);                                                                      
    bboxlo[0] = domain->sublo_lamda[0];
    bboxlo[1] = domain->sublo_lamda[1];
    bboxlo[2] = domain->sublo_lamda[2];
    bboxhi[0] = domain->subhi_lamda[0];
    bboxhi[1] = domain->subhi_lamda[1];
    bboxhi[2] = domain->subhi_lamda[2];
    bininvx *= domain->xprd;
    bininvy *= domain->yprd;
    bininvz *= domain->zprd;

See attached performs for examples/melt 500 atoms.

[stanmoore1]

@athomps is this for zero skew?

[stanmoore1]

This PR fixes the performance regression I saw for ReaxFF in PR #3740 on GPUs.

[sjplimp]

I don't think sorting in lamda coords is straightforward to do robustly. This is because the motion of atoms
and their closeness to each other is meausred in box coords. The point of sorting is to put atoms that are "close" in
box coords in the same (or nearby) bin so that when reordered they are also close in the X vector. "Closeness" in lamda
space is not the same as closeness in box space. In particular the distance metric between 2 particles is not
sqrt(dx^2 + dy^2 + dz^2) but is a function of the box skew, i.e. there is a different prefactor on dx vs dy vs dz.

[stanmoore1]

I agree with @sjplimp, I think we should just stick with converting the atom positions using lamda2x.

[athomps]

Yeah, I see. Atoms in the same skewed lamda bin can be arbitrarily far apart.

[athomps]

I approve.