Merge pull request #29 from chrisjonesBSU/add/unit_tests

Fix for `Lattice` box and `Tensile` simulation. Add more unit tests

hoomd_organics/base/system.py

@@ -520,7 +520,14 @@ def _build_system(self):
             layer.translate((self.x * i, 0, 0))
             system.add(layer)
         bounding_box = system.get_boundingbox()
-        x_len = bounding_box.lengths[0]
-        y_len = bounding_box.lengths[1]
+        # Add lattice constants to box lengths to account for boundaries
+        x_len = bounding_box.lengths[0] + self.x
+        y_len = bounding_box.lengths[1] + self.y
+        z_len = bounding_box.lengths[2] + 0.2
         self.set_target_box(x_constraint=x_len, y_constraint=y_len)
+        # Center the lattice in its box
+        system.box = mb.box.Box(np.array([x_len, y_len, z_len]))
+        system.translate_to(
+            (system.box.Lx / 2, system.box.Ly / 2, system.box.Lz / 2)
+        )
         return system

hoomd_organics/library/simulations/tensile.py

@@ -20,7 +20,7 @@ def __init__(
         gsd_write_freq=1e4,
         gsd_file_name="trajectory.gsd",
         log_write_freq=1e3,
-        log_file_name="sim_data.txt",
+        log_file_name="log.txt",
     ):
         super(Tensile, self).__init__(
             initial_state=initial_state,
@@ -34,16 +34,9 @@ def __init__(
             log_write_freq=log_write_freq,
             log_file_name=log_file_name,
         )
-        self.tensile_axis = tensile_axis.lower()
+        self.tensile_axis = np.asarray(tensile_axis)
         self.fix_ratio = fix_ratio
-        axis_array_dict = {
-            "x": np.array([1, 0, 0]),
-            "y": np.array([0, 1, 0]),
-            "z": np.array([0, 0, 1]),
-        }
-        axis_dict = {"x": 0, "y": 1, "z": 2}
-        self._axis_index = axis_dict[self.tensile_axis]
-        self._axis_array = axis_array_dict[self.tensile_axis]
+        self._axis_index = np.where(self.tensile_axis != 0)[0]
         self.initial_box = self.box_lengths_reduced
         self.initial_length = self.initial_box[self._axis_index]
         self.fix_length = self.initial_length * fix_ratio
@@ -88,7 +81,7 @@ def run_tensile(self, strain, kT, n_steps, period):
         )
         particle_puller = PullParticles(
             shift_by=shift_by / 2,
-            axis=self._axis_array,
+            axis=self.tensile_axis,
             neg_filter=self.fix_left,
             pos_filter=self.fix_right,
         )

hoomd_organics/tests/library/test_tensile.py

@@ -1,6 +1,29 @@
+import numpy as np
+
+from hoomd_organics.base.system import Lattice
+from hoomd_organics.library import OPLS_AA_PPS, PPS, Tensile
 from hoomd_organics.tests import BaseTest
 
 
 class TestTensileSimulation(BaseTest):
-    def test_tensile(self, polyethylene_system):
-        pass
+    def test_tensile(self):
+        pps = PPS(lengths=6, num_mols=32)
+        system = Lattice(
+            molecules=[pps],
+            force_field=[OPLS_AA_PPS()],
+            density=1.0,
+            r_cut=2.5,
+            x=1.2,
+            y=1.2,
+            n=4,
+            auto_scale=True,
+        )
+        tensile_sim = Tensile(
+            initial_state=system.hoomd_snapshot,
+            forcefield=system.hoomd_forcefield,
+            tensile_axis=(1, 0, 0),
+            log_write_freq=1e6,
+            gsd_write_freq=1e6,
+        )
+        tensile_sim.run_tensile(strain=0.05, kT=2.0, n_steps=1e3, period=10)
+        assert np.allclose(tensile_sim.strain, 0.05, 1e-4)

hoomd_organics/tests/modules/test_weld.py

@@ -1,8 +1,10 @@
 import os
 
 import gsd.hoomd
+import hoomd
 
 from hoomd_organics import Simulation
+from hoomd_organics.modules.utils import add_void_particles
 from hoomd_organics.modules.welding import Interface, SlabSimulation
 from hoomd_organics.tests.base_test import BaseTest
 
@@ -93,3 +95,22 @@ def test_weld_sim(self, polyethylene_system):
         )
         if os.path.isfile("restart.gsd"):
             os.remove("restart.gsd")
+
+    def test_void_particle(self, polyethylene_system):
+        init_snap = polyethylene_system.hoomd_snapshot
+        init_num_particles = init_snap.particles.N
+        init_types = init_snap.particles.types
+        void_snap, ff = add_void_particles(
+            init_snap,
+            polyethylene_system.hoomd_forcefield,
+            void_diameter=0.4,
+            num_voids=1,
+            void_axis=(1, 0, 0),
+            epsilon=1,
+            r_cut=0.4,
+        )
+        assert init_num_particles == void_snap.particles.N - 1
+        lj = [i for i in ff if isinstance(i, hoomd.md.pair.LJ)][0]
+        for p_type in init_types:
+            assert lj.params[(p_type, "VOID")]["sigma"] == 0.4
+            assert lj.params[(p_type, "VOID")]["epsilon"] == 1