unit tests; load only pdb or mol2

mbuild/compound.py

@@ -23,32 +23,21 @@
 
 
 def load(filename, relative_to_module=None, frame=-1, compound=None,
-         coords_only=False, **kwargs):
+         coords_only=False,**kwargs):
     """ """
-
     # Handle mbuild *.py files containing a class that wraps a structure file
     # in its own folder. E.g., you build a system from ~/foo.py and it imports
     # from ~/bar/baz.py where baz.py loads ~/bar/baz.pdb.
     if relative_to_module:
-        current_dir = os.path.dirname(os.path.realpath(
-            sys.modules[relative_to_module].__file__))
-        filename = os.path.join(current_dir, filename)
+        script_path = os.path.realpath(sys.modules[relative_to_module].__file__)
+        file_dir = os.path.dirname(script_path)
+        filename = os.path.join(file_dir, filename)
 
-    # This can return a md.Trajectory or a mb.Compound.
-    loaded = md.load(filename, **kwargs)
+    if compound is None:
+        compound = Compound()
 
-    if not compound:
-        if isinstance(loaded, Compound):
-            return loaded
-        else:
-            compound = Compound()
-
-    if isinstance(loaded, md.Trajectory):
-        compound.from_trajectory(loaded, frame=frame, coords_only=coords_only)
-    elif isinstance(loaded, Compound):  # Only updating coordinates.
-        assert compound.n_atoms == loaded.n_atoms
-        for atom, loaded_atom in zip(compound.atoms, loaded.atoms):
-            atom.pos = loaded_atom.pos
+    traj = md.load(filename, **kwargs)
+    compound.from_trajectory(traj, frame=frame, coords_only=coords_only)
     return compound
 
 
@@ -446,8 +435,7 @@ def from_trajectory(self, traj, frame=-1, coords_only=False):
         else:
             self.periodicity = np.array([0., 0., 0.])
 
-    def to_trajectory(self, show_ports=False, chain_types=None,
-                      residue_types=None, forcefield=None):
+    def to_trajectory(self, show_ports=False, chain_types=None, residue_types=None):
         """Convert to an md.Trajectory and flatten the compound.
 
         This also produces an object subclassed from MDTraj's Topology which
@@ -471,8 +459,7 @@ def to_trajectory(self, show_ports=False, chain_types=None,
         exclude = not show_ports
         atom_list = self.atom_list_by_name('*', exclude_ports=exclude)
 
-        top = self._to_topology(atom_list, chain_types=chain_types,
-                                residue_types=residue_types, forcefield=forcefield)
+        top = self._to_topology(atom_list, chain_types, residue_types)
 
         # Coordinates.
         xyz = np.ndarray(shape=(1, top.n_atoms, 3), dtype='float')
@@ -491,8 +478,7 @@ def to_trajectory(self, show_ports=False, chain_types=None,
         return md.Trajectory(xyz, top, unitcell_lengths=unitcell_lengths,
                              unitcell_angles=np.array([90, 90, 90]))
 
-    def _to_topology(self, atom_list, chain_types=None, residue_types=None,
-                     forcefield=None):
+    def _to_topology(self, atom_list, chain_types=None, residue_types=None):
         """Create a Topology from a Compound.
 
         Parameters
@@ -611,7 +597,7 @@ def save(self, filename, show_ports=False, forcefield=None, **kwargs):
                              'information are: {0}'.format(ff_formats))
 
         if saver:  # mBuild supported saver.
-            traj = self.to_trajectory(forcefield=forcefield, **kwargs)
+            traj = self.to_trajectory(**kwargs)
             return saver(filename, traj, show_ports=show_ports)
         else:  # MDTraj supported saver.
             traj = self.to_trajectory(show_ports=show_ports, **kwargs)

mbuild/examples/alkane/alkane.py

@@ -45,8 +45,8 @@ def __init__(self, n=3, cap_front=True, cap_end=True):
 
 def main():
     n = 5
-    alkane = Alkane(n=n, cap_front=True, cap_end=True)
-    alkane.visualize(show_ports=True)
+    return Alkane(n=n, cap_front=True, cap_end=True)
 
 if __name__ == "__main__":
-    main()
+    alkane = main()
+    alkane.visualize(show_ports=True)

mbuild/examples/alkane_monolayer/alkane_monolayer.py

@@ -29,7 +29,7 @@ def __init__(self, mask, tile_x=1, tile_y=1, chain_length=10):
 
         # Attach chains to specified binding sites. Other sites get a hydrogen.
         mb.apply_mask(host=self.tiled_surface, guest=alkylsilane, mask=mask,
-                   backfill=hydrogen)
+                      backfill=hydrogen)
 
 
 def main():

mbuild/examples/bilayer/bilayer.py

@@ -225,9 +225,10 @@ def main():
                       solvent_per_lipid=20, mirror=False)
 
     bilayer.save(filename='bilayer.pdb')
-
-    #import os
-    #os.system('vmd -e vis.vmd')
+    return bilayer
 
 if __name__ == "__main__":
     main()
+    import os
+    os.system('vmd -e vis.vmd')
+

mbuild/examples/ethane/ethane.py

@@ -21,4 +21,5 @@ def main():
 
 if __name__ == "__main__":
     ethane = main()
+    import ipdb; ipdb.set_trace()
     ethane.visualize(show_ports=True)

mbuild/examples/methane/methane.py

@@ -26,3 +26,4 @@ def main():
 
 if __name__ == "__main__":
     main()
+    methane.visualize()

mbuild/packing.py

@@ -55,7 +55,8 @@ def fill_box(compound, n_compounds, box, overlap=0.2):
 
     compound_pdb = tempfile.mkstemp(suffix='.pdb')[1]
     compound.save(compound_pdb)
-    filled_pdb = tempfile.mkstemp(suffix='.pdb')[1]
+    #filled_pdb = tempfile.mkstemp(suffix='.pdb')[1]
+    filled_pdb = 'filled.pdb'
 
     # In angstroms for packmol.
     box_lengths = box.lengths * 10
@@ -67,6 +68,8 @@ def fill_box(compound, n_compounds, box, overlap=0.2):
 
     proc = Popen(PACKMOL, stdin=PIPE, stdout=PIPE, stderr=PIPE, universal_newlines=True)
     out, err = proc.communicate(input=input_text)
+    if err:
+        _packmol_error(out, err)
 
     # Create the topology and update the coordinates.
     filled = Compound()
@@ -95,11 +98,15 @@ def solvate(solute, solvent, n_solvent, box, overlap=0.2):
     if not PACKMOL:
         raise IOError("Packmol not found")
 
+    if isinstance(box, (list, tuple)):
+        box = Box(lengths=box)
+
     solute_pdb = tempfile.mkstemp(suffix='.pdb')[1]
     solute.save(solute_pdb)
     solvent_pdb = tempfile.mkstemp(suffix='.pdb')[1]
     solvent.save(solvent_pdb)
-    solvated_pdb = tempfile.mkstemp(suffix='.pdb')[1]
+    #solvated_pdb = tempfile.mkstemp(suffix='.pdb')[1]
+    solvated_pdb = 'solvated.pdb'
 
     # In angstroms for packmol.
     box_lengths = box.lengths * 10
@@ -113,6 +120,8 @@ def solvate(solute, solvent, n_solvent, box, overlap=0.2):
 
     proc = Popen(PACKMOL, stdin=PIPE, stdout=PIPE, stderr=PIPE, universal_newlines=True)
     out, err = proc.communicate(input=input_text)
+    if err:
+        _packmol_error(out, err)
 
     # Create the topology and update the coordinates.
     solvated = Compound()
@@ -121,3 +130,9 @@ def solvate(solute, solvent, n_solvent, box, overlap=0.2):
         solvated.add(deepcopy(solvent))
     solvated.update_coordinates(solvated_pdb)
     return solvated
+
+def _packmol_error(out, err):
+    with open('log.txt', 'w') as log_file, open('err.txt', 'w') as err_file:
+        log_file.write(out)
+        err_file.write(err)
+    raise RuntimeError("PACKMOL failed. See 'err.txt' and 'log.txt'")

mbuild/port.py

@@ -4,6 +4,7 @@
 
 from mbuild.atom import Atom
 from mbuild.compound import Compound
+from mbuild.coordinate_transform import rotate_around_z
 
 
 class Port(Compound):
@@ -38,7 +39,6 @@ def __init__(self, anchor=None):
 
         down = deepcopy(up)
 
-        from mbuild.coordinate_transform import rotate_around_z
         rotate_around_z(down, np.pi)
 
         self.add(up, 'up')

mbuild/tests/base_test.py

@@ -30,3 +30,13 @@ def methane(self):
     def h2o(self):
         from mbuild.components.small_groups.h2o import H2O
         return H2O()
+
+    @pytest.fixture
+    def ch2(self):
+        from mbuild.components.small_groups.ch2 import CH2
+        return CH2()
+
+    @pytest.fixture
+    def betacristobalite(self):
+        from mbuild.components.surfaces.betacristobalite import Betacristobalite
+        return Betacristobalite()

mbuild/tests/test_mask.py

@@ -0,0 +1,35 @@
+import pytest
+
+import mbuild as mb
+from mbuild.tests.base_test import BaseTest
+
+
+class TestMask(BaseTest):
+
+    @pytest.mark.skipif(True, reason='Needs implementing!')
+    def test_apply_mask(self):
+        pass
+
+    def test_random_2d(self):
+        mask = mb.random_mask_2d(100)
+        assert len(mask) == 100
+
+    def test_random_3d(self):
+        mask = mb.random_mask_3d(100)
+        assert len(mask) == 100
+
+    def test_grid_2d(self):
+        mask = mb.grid_mask_2d(10, 5)
+        assert len(mask) == 50
+
+    def test_grid_3d(self):
+        mask = mb.grid_mask_3d(10, 5, 2)
+        assert len(mask) == 100
+
+    def test_sphere(self):
+        mask = mb.sphere_mask(100)
+        assert len(mask) == 100
+
+    def test_disk(self):
+        mask = mb.disk_mask(100)
+        assert len(mask) == 100

mbuild/tests/test_packing.py

@@ -0,0 +1,15 @@
+import mbuild as mb
+from mbuild.tests.base_test import BaseTest
+
+
+class TestPacking(BaseTest):
+
+    def test_fill_box(self, ethane):
+        filled = mb.fill_box(ethane, n_compounds=20, box=[2, 2, 2])
+        assert filled.n_atoms == 20 * 8
+        assert filled.n_bonds == 20 * 7
+
+    def test_solvate(self, ethane, h2o):
+        solvated = mb.solvate(ethane, h2o, n_solvent=100, box=[4, 4, 4])
+        assert solvated.n_atoms == 8 + 100 * 3
+        assert solvated.n_bonds == 7 + 100 * 2

mbuild/tests/test_polymer.py

@@ -0,0 +1,11 @@
+import mbuild as mb
+from mbuild.tests.base_test import BaseTest
+
+
+class TestPolymer(BaseTest):
+
+    def test_polymer(self, ch2):
+        n = 6
+        c6 = mb.Polymer(ch2, n=n)
+        assert c6.n_atoms == n * 3
+        assert c6.n_bonds == n * 2 + (n - 1)

mbuild/tests/test_tiled_compound.py

@@ -0,0 +1,31 @@
+import pytest
+
+import mbuild as mb
+from mbuild.tests.base_test import BaseTest
+
+
+class TestTiledCompound(BaseTest):
+
+    def test_2d_replication(self, betacristobalite):
+        nx = 2
+        ny = 3
+        nz = 1
+        tiled = mb.TiledCompound(betacristobalite, [nx, ny, nz])
+        assert tiled.n_atoms == 1800 * nx * ny * nz
+        assert tiled.n_bonds == 2300 * nx * ny * nz
+
+    def test_incorrect_periodicity(self, betacristobalite):
+        nx = 2
+        ny = 3
+        nz = 2
+        with pytest.raises(ValueError):
+            mb.TiledCompound(betacristobalite, [nx, ny, nz])
+
+    def test_negative_periodicity(self, betacristobalite):
+        nx = -2
+        ny = 3
+        nz = 2
+        with pytest.raises(ValueError):
+            mb.TiledCompound(betacristobalite, [nx, ny, nz])
+
+

mbuild/tiled_compound.py

@@ -33,15 +33,15 @@ def __init__(self, tile, n_tiles=None, kind=None):
         if not n_tiles:
             n_tiles = (1, 1, 1)
         n_x, n_y, n_z = n_tiles
-        assert n_x > 0 and n_y > 0 and n_z > 0, "Number of tiles must be positive."
+        if not n_x > 0 and n_y > 0 and n_z > 0:
+            raise ValueError("Number of tiles must be positive.")
 
         # Check that the tile is periodic in the requested dimensions.
-        if n_x != 1:
-            assert tile.periodicity[0] != 0, "Tile not periodic in x dimension."
-        if n_y != 1:
-            assert tile.periodicity[1] != 0, "Tile not periodic in y dimension."
-        if n_z != 1:
-            assert tile.periodicity[2] != 0, "Tile not periodic in z dimension."
+        if ((n_x != 1 and tile.periodicity[0] == 0) or
+                (n_y != 1 and tile.periodicity[1] == 0) or
+                (n_z != 1 and tile.periodicity[2] == 0)):
+            raise ValueError("Tile not periodic in at least one of the "
+                             "specified dimensions.")
 
         self.tile = tile
         self.n_x = n_x

mbuild/utils/io.py

@@ -1,10 +1,6 @@
-from __future__ import print_function
-
 import os
 from pkg_resources import resource_filename
 
-from mbuild.compound import load
-
 
 def get_fn(name):
     """Get the full path to one of the reference files shipped for utils