Support for periodic boundary condition!

[mushroomfire]

Can you privode support of periodic boundary condition when calculating distances? Just like the boxsize parameter in scipy.spatial.Kdtree. Thanks a lot.

[djhoese]

I'll be honest I'm not exactly sure how this parameter works and I also don't know the low-level details of the pykdtree algorithm. The maintainers of this library don't really have time for new features like this, but if someone can read the source code, figure out how it can be added, and make a pull request, we will happily review it and merge it (as long as it doesn't hurt existing performance).

[ivan-pi]

The idea is just to restrict the particle coordinates (see https://en.wikipedia.org/wiki/Periodic_boundary_conditions).

When calculating distances, you need to modify the calculation as follows:

if (periodic_x) then
  dx = x(j) - x(i)
  if (dx >   x_size * 0.5) dx = dx - x_size
  if (dx <= -x_size * 0.5) dx = dx + x_size
end if

(assuming the periodic box is centered at the coordinate system origin)

[djhoese]

Thanks for the wikipedia reference and the example code. Modifying the .mako template with these changes and then making a pull request would be the next best step.

[mushroomfire]

Hello, @djhoese Actually, we can make it more generalizable for both rectangular and triclinic box. Let's consider the box below, where each row represents a box vector:

box = np.array([[a1, a2, a3], [b1, b2, b3], [c1, c2, c3]])

The rectangular box is like below:

box = np.array([[a, 0, 0], [0, b, 0], [0, 0, c]])

The original boxsize is the length of the rectangular box boxsize=[a, b, c].

For two kinds of box, a distance vector rij = np.array([x, y, z]) can be wrapped as below:

def pbc(rij, box):
      nz = rij[2] / box[2][2]
      ny = (rij[1] - nz * box[2][1]) / box[1][1]
      nx = (rij[0] - ny * box[1][0] - nz * box[2][0]) / box[0][0]
      n =[nx, ny, nz]
      for i in range(3):
            if n[i] > 0.5:
                n[i] -= 1
            elif n[i] < -0.5:
                n[i] += 1
      return n[0] * box[0] + n[1] * box[1] + n[2] * box[2]

Maybe change the distance equation in codes can achieve this?

Thank you very much.

[djhoese]

I'm not sure I understand, but I'm also not sure I need to. You'll be able to get more feedback by creating a pull request so that the other maintainers who have even more limited time for this project than I do will be able to see concrete examples of what needs to change, what the benefits are, and can provide feedback.