Add basic method to set up Monolayer and Solvated Monolayer (#2)

* add basic method to set up monolayer and solvated monolayer

* start adding simple tests

* Update docstring for most calsses

environment.yml

@@ -11,4 +11,5 @@ dependencies:
   - py3dmol
   - python=3.8
   - rdkit
-  - scipy
+  - scipy
+  - pytest

surface_coatings/__init__.py

@@ -1 +1,2 @@
-from .monolayer import Monolayer
+from .monolayer import Monolayer, DualMonolayer
+from .solvated_monolayer import SolvatedMonolayer, SolvatedDualMonolayer

surface_coatings/chains/alkylsilane.py

@@ -1,3 +1,4 @@
+"""Routine to create alkylsilane chain."""
 import mbuild as mb
 from mbuild.lib.recipes import Alkane
 from mbuild.lib.moieties import Silane
@@ -35,4 +36,4 @@ def __init__(self, chain_length=17, terminal_group="methyl"):
         self.add(silane, 'silane')
         mb.force_overlap(self['silane'], self['silane']['up'], self['alkane']['down'])
 
-        self.add(silane['down'], 'down', containment=False)
+        self.add(silane['down'], 'down', con tainment=False)

surface_coatings/chains/vbc_polymer.py

@@ -1,3 +1,4 @@
+"""Routine to create VBC polymer."""
 import mbuild as mb
 from mbuild.lib.recipes import Polymer
 import numpy as np
@@ -9,22 +10,36 @@
 class VBCPolymer(mb.Compound):
     def __init__(self, monomers=[Methacrylate(), SBMA(), AzPMA(), TriazoleBiotin()], n=1,
                  sequence='AABCBBD', port_labels=('up', 'down')):
+        """This is a general purpose to create a Vinylbenzyl Polymer
+
+        Parameters
+        ----------
+        monomers: list of monomer
+            List of monomers to be connected together
+        n: int, optional, default=1
+            Number of repeat for the polymer
+        sequence: str
+            The sequence of all the monomer, corresponding to the monomers list
+        port_labels: tuple, optional, default=('up', 'down')
+            The list ports of the monomers. The VBC will connect to the first port,
+            and the last port will be capped by a Hydrogen
+        """
         super(VBCPolymer, self).__init__()
         copolymer = Polymer(monomers=monomers)
         copolymer.build(n=n, sequence=sequence, add_hydrogens=False)
-        self.add(copolymer)
+        self.add(copolymer, label="Polymer")
         vbc = MethylStyrene()
-        self.add(vbc)
+        self.add(vbc, label="VBC")
         mb.force_overlap(vbc,
                          vbc['up'],
-                         copolymer['up'])
+                         copolymer[port_labels[0]])
         for port in self.all_ports():
             if port.access_labels:
                 self.labels['up'] = port
             else:
                 self.labels['down'] = port
 
-        orientation = self['Polymer[0]']['down'].direction
+        orientation = self['Polymer']['down'].direction
         for port in self.available_ports():
             if not all(np.isclose(port.direction, orientation)):
                 port.update_orientation(-orientation)

surface_coatings/monolayer.py

@@ -1,6 +1,4 @@
-from copy import deepcopy
-from warnings import warn
-
+"""Routines to create (dual) monolayer systems."""
 import numpy as np
 
 import mbuild as mb
@@ -9,11 +7,24 @@
 
 
 class Monolayer(mb.Compound):
-    """A general monolayer recipe.
+    """A surface coated by a monolayer.
 
     Parameters
     ----------
-
+    surface: mb.Compound
+        The surface with ports at its surface.
+    chain: mb.Compound
+        The chain that to be attached to the surfaces.
+    n_chains: int
+        The number of chains to be attached.
+    backfill: mb.Compound, optional, default=H()
+        Compound used to backfill leftover ports (after all chains have been attached.
+    tile_x, tile_y: int, optional, default= 1, 1
+        The number of surface tiles.
+    rotate: bool, optional, default=True
+        Options to rotate the chain randomly.
+    seed: int, optional, default= 12345
+        Random seed used for any subprocess.
     """
     def __init__(self, surface, chain, n_chains, backfill=H(), tile_x=1, tile_y=1, rotate=True, seed=12345, **kwargs):
         super(Monolayer, self).__init__()
@@ -37,4 +48,29 @@ def __init__(self, surface, chain, n_chains, backfill=H(), tile_x=1, tile_y=1, r
                 rotation = np.random.random() * np.pi * 2.0
                 chain.spin(rotation, [0, 0, 1])
 
-        self.periodicity = [True, True, False]
+        self.periodicity = surface.periodicity
+
+
+class DualMonolayer(mb.Compound):
+    """A dual-monolayer system.
+
+    Parameters
+    ----------
+    top: mb.Compound
+        The top surfaces of the dual-monolayer.
+    bottom: mb.Compound
+        The bottom surface of the dual-monolayer.
+    separation: float, optional, default=0.8
+        The separation between the two surfaces.
+    """
+    def __init__(self, top, bottom, separation=0.8):
+        super(DualMonolayer, self).__init__()
+        top.spin(np.pi, around=[0, 1, 0])
+        top_box = top.get_boundingbox()
+        bot_box = bottom.get_boundingbox()
+
+        z_val = bot_box.lengths[2]
+
+        top.translate([0, 0, z_val + separation])
+        self.add(top, label="top_monolayer")
+        self.add(bottom, label="bottom_monolayer")

surface_coatings/solvated_monolayer.py

@@ -0,0 +1,68 @@
+"""Routines to construct solvated (dual) monolayer."""
+import numpy as np
+import mbuild as mb
+from mbuild.lib.moieties import H2O
+
+
+class SolvatedMonolayer(mb.Compound):
+    """Solvated monolayer system.
+
+    Parameters
+    ----------
+    monolayer: mb.Compound
+        The monolayer to be solvated.
+    solvent: mb.Compound, optional, default=H2O()
+        The solvent to be used.
+    n_solvents: int, optional, default=1000
+        The number of solvent compounds used.
+    solvent_box_height: float, optional, default=5
+        The height of the solvent box. The base of the box adapt those of the monolayer/surface.
+    seed: int, optional, default=12345
+        Random seed used for any subprocess.
+    """
+    def __init__(self, monolayer, solvent=H2O(), n_solvents=1000, solvent_box_height=5, seed=12345):
+        super(SolvatedMonolayer, self).__init__()
+        monolayer_box_lengths = monolayer.get_boundingbox().lengths
+        solvent_box = [monolayer_box_lengths[0],
+                       monolayer_box_lengths[1],
+                       solvent_box_height]
+        box_of_solvent = mb.fill_box(compound=solvent,
+                                     box=solvent_box,
+                                     n_compounds=n_solvents)
+        box_of_solvent.translate([0, 0, monolayer_box_lengths[2]])
+        self.add(monolayer, label="monolayer")
+        self.add(box_of_solvent, label="solvent")
+
+
+class SolvatedDualMonolayer(mb.Compound):
+    """Solvated dual-monolayer system.
+
+    Parameters
+    ----------
+    dual_monolayer: mb.Compound
+        The dual-monolayer system to be solvated.
+    solvent: mb.Compound, optional, default=H2O()
+        The solvent compound to be used.
+    n_solvents: int, optional, n=1000
+        The number of solvent molecules to be used.
+    seed: int, optional, default=12345
+        Random seed used for any subprocess.
+    """
+    def __init__(self, dual_monolayer, solvent=H2O(), n_solvents=1000, seed=12345):
+        super(SolvatedDualMonolayer, self).__init__()
+        top_monolayer = dual_monolayer["top_monolayer"]
+        top_monolayer_box_lengths = top_monolayer.get_boundingbox().lengths
+        bottom_monolayer = dual_monolayer["bottom_monolayer"]
+        bottom_monolayer_box_lengths = bottom_monolayer.get_boundingbox().lengths
+        dual_monolayer_box_lengths = dual_monolayer.get_boundingbox().lengths
+
+        separation = dual_monolayer_box_lengths[2] - (top_monolayer_box_lengths[2] + bottom_monolayer_box_lengths[2])
+        solvent_box = [bottom_monolayer_box_lengths[0],
+                       bottom_monolayer_box_lengths[1],
+                       separation]
+        box_of_solvent = mb.fill_box(compound=solvent,
+                                     box=solvent_box,
+                                     n_compounds=n_solvents)
+        box_of_solvent.translate([0, 0, bottom_monolayer_box_lengths[2]])
+        self.add(dual_monolayer, label="monolayers")
+        self.add(box_of_solvent, label="solvents")

surface_coatings/surfaces/__init__.py

@@ -0,0 +1,2 @@
+from .silicon_interface import SiliconInterface
+from .silica_interface_carve import SilicaInterfaceCarve

surface_coatings/surfaces/silica_interface_carve.py

@@ -1,3 +1,4 @@
+"""Routine to create carve amorphous silica interface."""
 from __future__ import division
 
 import math
@@ -167,8 +168,3 @@ def _adjust_stoichiometry(self):
             O1 = random.choice(bottom_Os)
             bottom_Os.remove(O1)
             self.remove(O1)
-
-if __name__ == "__main__":
-    from mbuild.lib.bulk_materials import AmorphousSilica
-    silica_interface = mb.SilicaInterfaceCarve(bulk_silica=AmorphousSilica(), thickness=1.2)
-    silica_interface.save('silica_interface.mol2', show_ports=True)

surface_coatings/surfaces/silicon_interface.py

@@ -1,9 +1,16 @@
+"""Routine to create silica surface."""
 import numpy as np
 import mbuild as mb
-from surface_coatings.monomers.methylstyrene import MethylStyrene
 
 
 class CrystalineSilicon(mb.Compound):
+    """A block of crystaline silicon.
+
+    Parameters
+    ----------
+    x, y, z: float, default= 10, 10, 10
+        Dimension of the silicon block.
+    """
     def __init__(self, x=10, y=10, z=10):
         super(CrystalineSilicon, self).__init__()
         # define all necessary lattice parameters
@@ -30,6 +37,17 @@ def __init__(self, x=10, y=10, z=10):
 
 
 class SiliconInterface(mb.Compound):
+    """A surface made up of crystaline silicon.
+
+    Parameters
+    ----------
+    silicon: CrystalineSilicon
+        The crystaline silicon building block.
+    tile_x, tile_y: int
+        The number of tiles to build out the surface.
+    seed: int, optional, default=12345
+        Random seed used in some subprocess.
+    """
     def __init__(self, silicon=CrystalineSilicon(x=10, y=10, z=2), tile_x=1, tile_y=1, seed=12345):
         super(SiliconInterface, self).__init__()
         tiled_compound = mb.lib.recipes.TiledCompound(mb.clone(silicon), n_tiles=(tile_x, tile_y, 1))
@@ -39,6 +57,7 @@ def __init__(self, silicon=CrystalineSilicon(x=10, y=10, z=2), tile_x=1, tile_y=
         self.periodicity = silicon.periodicity
 
     def _identify_surface_sites(self):
+        """Method to identify and add port to silicon at the surface."""
         for particle in list(self.particles()):
             if np.isclose(particle.pos[2], 0):
                 label = f'Si_{len(self.referenced_ports())}'

surface_coatings/tests/test_surfaces.py

@@ -0,0 +1,12 @@
+import mbuild as mb
+from surface_coatings.surfaces import SiliconInterface, SilicaInterfaceCarve
+
+
+class TestSystem(object):
+    def test_silicon_interface(self):
+        silicon_surface = SiliconInterface()
+        assert silicon_surface.periodicity == (True, True, False)
+
+    def test_silica_interface_carve(self):
+        silica_surface = SilicaInterfaceCarve)
+        assert silica_surface.periodicity == (True, True, False)