freud.locality.PeriodicBuffer: Unit Cell RDF

The PeriodicBuffer class is meant to replicate points beyond a single image while respecting box periodicity. This example demonstrates how we can use this to compute the radial distribution function from a sample crystal’s unit cell.

%matplotlib inline
import freud
import matplotlib.pyplot as plt
import numpy as np

Here, we create a box to represent the unit cell and put two points inside. We plot the box and points below.

box = freud.box.Box(Lx=2, Ly=2, xy=np.sqrt(1 / 3), is2D=True)
points = np.array([[-0.5, -0.5, 0], [0.5, 0.5, 0]])
system = freud.AABBQuery(box, points)
system.plot(ax=plt.gca())
plt.show()

Next, we create a PeriodicBuffer instance and have it compute the “buffer” points that lie outside the first periodicity. These positions are stored in the buffer_points attribute. The corresponding buffer_ids array gives a mapping from the index of the buffer particle to the index of the particle it was replicated from, in the original array of points. Finally, the buffer_box attribute returns a larger box, expanded from the original box to contain the replicated points.

pbuff = freud.locality.PeriodicBuffer()
pbuff.compute(system=(box, points), buffer=6, images=True)
print(pbuff.buffer_points[:10], "...")

[[ 0.65470022  1.5         0.        ]
 [ 1.80940032  3.5         0.        ]
 [ 2.96410179  5.5         0.        ]
 [-3.96410131 -6.5         0.        ]
 [-2.80940104 -4.49999952  0.        ]
 [-1.65470016 -2.50000048  0.        ]
 [ 1.50000024 -0.5         0.        ]
 [ 2.65470076  1.5         0.        ]
 [ 3.80940032  3.5         0.        ]
 [ 4.96410179  5.5         0.        ]] ...

Below, we plot the original unit cell and the replicated buffer points and buffer box.

system.plot(ax=plt.gca())
plt.scatter(pbuff.buffer_points[:, 0], pbuff.buffer_points[:, 1])
pbuff.buffer_box.plot(ax=plt.gca(), linestyle="dashed", color="gray")
plt.show()

Finally, we can plot the radial distribution function (RDF) of this replicated system, using a value of r_max that is larger than the size of the original box. This allows us to see the interaction of the original points with their replicated neighbors from the buffer.

rdf = freud.density.RDF(bins=250, r_max=5)
rdf.compute(system=(pbuff.buffer_box, pbuff.buffer_points), query_points=points)
rdf.plot(ax=plt.gca())
plt.show()