From 28c645b905c004199021d8e55ff255cf535acc36 Mon Sep 17 00:00:00 2001 From: jac16 Date: Tue, 28 Nov 2023 09:48:07 -0500 Subject: [PATCH 01/11] add make_periodic option to mbuild --- mbuild/lib/recipes/polymer.py | 57 +++++++++++++++++++++++------------ 1 file changed, 37 insertions(+), 20 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 3b1eb623b..2ed1b4bba 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -107,7 +107,7 @@ def end_groups(self): """ return self._end_groups - def build(self, n, sequence="A", add_hydrogens=True): + def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): """Connect one or more components in a specified sequence. Uses the compounds that are stored in Polymer.monomers and @@ -126,11 +126,18 @@ def build(self, n, sequence="A", add_hydrogens=True): For example, 'AB' where 'A'corresponds to the first compound added to Polymer.monomers and 'B' to the second compound. add_hydrogens : bool, default True - If True and an end group compound is None, then the head or tail + If True and an ``end_groups`` compound is None, then the head or tail of the polymer will be capped off with hydrogen atoms. If end group compounds exist, then they will be used. If False and an end group compound is None, then the head or tail port will be exposed in the polymer. + make_periodic : bool, default False + If True and an ``end_groups`` compound is None, then the head and tail + will be forced into an overlap with a periodicity along the z-axis. + If end group compounds exist, then there will be a warning. However + ``add_hydrogens`` will simply be overwritten. + If False, ``end_groups`` compound is None, and ``add_hydrogens`` is + False then the head or tail port will be exposed in the polymer. """ if n < 1: raise ValueError("n must be 1 or more") @@ -182,25 +189,35 @@ def build(self, n, sequence="A", add_hydrogens=True): head_tail = [head_port, tail_port] - for i, compound in enumerate(self._end_groups): - if compound is not None: - if self._headtail[i] is not None: - head_tail[i].update_separation(self._headtail[i]) - self.add(compound) - force_overlap(compound, compound.labels["up"], head_tail[i]) - else: - if add_hydrogens: - hydrogen = H() - # Defaut to 1/2 H-C bond len - head_tail[i].update_separation(0.0547) - hydrogen["up"].update_separation(0.0547) - self.add(hydrogen) - force_overlap(hydrogen, hydrogen["up"], head_tail[i]) + if make_periodic: + if all(self._end_groups == None): + raise ValueError( + "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." + ) + self.periodic = [False, False, True] + force_overlap(self, head_tail[0], head_tail[1]) + else: + for i, compound in enumerate(self._end_groups): + if compound is not None: + if self._headtail[i] is not None: + head_tail[i].update_separation(self._headtail[i]) + self.add(compound) + force_overlap(compound, compound.labels["up"], head_tail[i]) else: - # if None, hoist port to polymer level - self.add( - head_tail[i], self._port_labels[i], containment=False - ) + if add_hydrogens: + hydrogen = H() + # Defaut to 1/2 H-C bond len + head_tail[i].update_separation(0.0547) + hydrogen["up"].update_separation(0.0547) + self.add(hydrogen) + force_overlap(hydrogen, hydrogen["up"], head_tail[i]) + else: + # if None, hoist port to polymer level + self.add( + head_tail[i], + self._port_labels[i], + containment=False, + ) for port in self.all_ports(): if port not in self.available_ports(): From 7009c98a270fbd7afb1638f9481f11c7d17da9d8 Mon Sep 17 00:00:00 2001 From: jac16 Date: Tue, 28 Nov 2023 13:03:43 -0500 Subject: [PATCH 02/11] Make periodic polymer --- mbuild/lib/recipes/polymer.py | 26 ++++++++++++++++++++++++-- 1 file changed, 24 insertions(+), 2 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 2ed1b4bba..042213325 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -190,12 +190,34 @@ def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): head_tail = [head_port, tail_port] if make_periodic: - if all(self._end_groups == None): + if np.any([x is not None for x in self._end_groups]): raise ValueError( "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." ) - self.periodic = [False, False, True] + self.translate(-head_port.pos) + vector0 = np.array([0, 0, 1]) + vector = tail_port.pos - head_port.pos + vector /= np.linalg.norm(vector) + dot = np.dot(vector0, vector) + if dot > 0: + theta = np.arccos(dot) + else: + theta = -(np.pi - np.arccos(dot)) + self.spin(theta, np.cross(vector0, vector), anchor=head_port) + self.periodic = list( + ~np.isclose( + head_port.pos - tail_port.pos, 0.0, np.finfo(float).eps + ) + ) force_overlap(self, head_tail[0], head_tail[1]) + if not np.all( + np.isclose(self.periodic, [0, 0, 1], np.finfo(float).eps) + ): + raise ValueError( + "Periodic polymer could not be aligned: {} {}".format( + head_port.pos, tail_port.pos + ) + ) else: for i, compound in enumerate(self._end_groups): if compound is not None: From c281341883773f53e0b6f6255a1a7f8b7473104b Mon Sep 17 00:00:00 2001 From: jac16 Date: Fri, 8 Mar 2024 15:50:53 -0500 Subject: [PATCH 03/11] Add test for periodic polymer --- mbuild/lib/recipes/polymer.py | 7 ++++++- mbuild/tests/test_monolayer.py | 29 +++++++++++++++++++++++++++++ 2 files changed, 35 insertions(+), 1 deletion(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 042213325..904db8ba7 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -83,7 +83,12 @@ def __init__(self, monomers=None, end_groups=None): super(Polymer, self).__init__() self._monomers = monomers or [] self._end_groups = end_groups or [None, None] - if len(self._end_groups) != 2: + if not isinstance(self._end_groups, list): + raise ValueError( + "Please provide two end groups in a list; " + f"you provided {self._end_groups}" + ) + elif len(self._end_groups) != 2: raise ValueError( "Please provide two end groups; " f"you provided {len(self._end_groups)}" diff --git a/mbuild/tests/test_monolayer.py b/mbuild/tests/test_monolayer.py index b1460aacf..214e37f85 100644 --- a/mbuild/tests/test_monolayer.py +++ b/mbuild/tests/test_monolayer.py @@ -1,3 +1,5 @@ +import pytest + import mbuild as mb from mbuild.lib.atoms import H from mbuild.lib.recipes import Monolayer, Polymer @@ -45,6 +47,33 @@ def test_pattern_kwargs(self, ch2): assert monolayer.n_particles == 2000 + chains * 29 assert monolayer.n_bonds == 2500 + chains * 29 + def test_periodic_pattern(self, ch2): + # Make Periodic without Hydrogen Conflict + chain = mb.recipes.Polymer(monomers=[ch2]) + chain.build(n=10, add_hydrogens=False, make_periodic=True) + assert not chain.all_ports() + + bonded_atoms = [ + x.name for x in list(chain["monomer[0]"][0].direct_bonds()) + ] + assert bonded_atoms.count("H") == 2 + assert bonded_atoms.count("C") == 2 + + # Make Periodic with Hydrogen Conflict + chain2 = mb.recipes.Polymer(monomers=[ch2]) + chain2.build(n=10, add_hydrogens=True, make_periodic=True) + assert not chain2.all_ports() + + bonded_atoms = [ + x.name for x in list(chain2["monomer[0]"][0].direct_bonds()) + ] + assert bonded_atoms.count("H") == 2 + assert bonded_atoms.count("C") == 2 + + with pytest.raises(ValueError): + chain3 = mb.recipes.Polymer(monomers=[ch2], end_groups=ch2) + chain3.build(n=10, add_hydrogens=True, make_periodic=True) + def test_mixed_monolayer(self, ch2): n = 8 m = 8 From 881f03a53d24cfb0263b24ae82c2df45dd913857 Mon Sep 17 00:00:00 2001 From: jac16 Date: Mon, 11 Mar 2024 14:50:45 -0400 Subject: [PATCH 04/11] Use bond length as a metric of periodicity --- mbuild/lib/recipes/polymer.py | 39 +++++++++++++++++++---------------- 1 file changed, 21 insertions(+), 18 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 904db8ba7..f6c81155b 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -6,7 +6,7 @@ from mbuild import clone from mbuild.compound import Compound -from mbuild.coordinate_transform import force_overlap +from mbuild.coordinate_transform import force_overlap, z_axis_transform from mbuild.lib.atoms import H from mbuild.port import Port from mbuild.utils.validation import assert_port_exists @@ -199,28 +199,31 @@ def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): raise ValueError( "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." ) - self.translate(-head_port.pos) - vector0 = np.array([0, 0, 1]) - vector = tail_port.pos - head_port.pos - vector /= np.linalg.norm(vector) - dot = np.dot(vector0, vector) - if dot > 0: - theta = np.arccos(dot) - else: - theta = -(np.pi - np.arccos(dot)) - self.spin(theta, np.cross(vector0, vector), anchor=head_port) - self.periodic = list( - ~np.isclose( - head_port.pos - tail_port.pos, 0.0, np.finfo(float).eps - ) + z_axis_transform( + compound=self, + point_on_z_axis=head_port.pos, + point_on_zx_plane=tail_port.pos, + ) + bond_length = 2 * np.max( + [ + np.sqrt(np.sum(np.square(atm1.pos - atm2.pos))) + for (atm1, atm2) in self.root.bond_graph.edges( + self["monomer[0]"] + ) + ] ) - force_overlap(self, head_tail[0], head_tail[1]) + self.periodic = [ + x > bond_length for x in head_port.pos - tail_port.pos + ] + force_overlap(self, head_port, tail_port) if not np.all( np.isclose(self.periodic, [0, 0, 1], np.finfo(float).eps) ): raise ValueError( - "Periodic polymer could not be aligned: {} {}".format( - head_port.pos, tail_port.pos + "Periodic polymer could not be aligned:\n Head Coord {}\n Tail Coord {}\n Difference {}".format( + head_port.pos, + tail_port.pos, + tail_port.pos - head_port.pos, ) ) else: From fda155a6348a7b09b55dc18ab8c82af321945999 Mon Sep 17 00:00:00 2001 From: Jennifer A Clark Date: Thu, 14 Mar 2024 18:48:46 -0400 Subject: [PATCH 05/11] Update mbuild/lib/recipes/polymer.py Co-authored-by: CalCraven <54594941+CalCraven@users.noreply.github.com> --- mbuild/lib/recipes/polymer.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index f6c81155b..eaf4ca807 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -206,7 +206,7 @@ def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): ) bond_length = 2 * np.max( [ - np.sqrt(np.sum(np.square(atm1.pos - atm2.pos))) + np.linalg.norm(atm1.pos - atm2.pos) for (atm1, atm2) in self.root.bond_graph.edges( self["monomer[0]"] ) From 358104868c6daade8390a439fc207c78f692572a Mon Sep 17 00:00:00 2001 From: Jennifer A Clark Date: Thu, 14 Mar 2024 18:48:46 -0400 Subject: [PATCH 06/11] Separate periodic connection into a method --- mbuild/lib/recipes/polymer.py | 105 +++++++++++++++++++++------------- 1 file changed, 65 insertions(+), 40 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index f6c81155b..3918792a9 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -6,7 +6,12 @@ from mbuild import clone from mbuild.compound import Compound -from mbuild.coordinate_transform import force_overlap, z_axis_transform +from mbuild.coordinate_transform import ( + force_overlap, + x_axis_transform, + y_axis_transform, + z_axis_transform, +) from mbuild.lib.atoms import H from mbuild.port import Port from mbuild.utils.validation import assert_port_exists @@ -112,7 +117,7 @@ def end_groups(self): """ return self._end_groups - def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): + def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): """Connect one or more components in a specified sequence. Uses the compounds that are stored in Polymer.monomers and @@ -136,12 +141,14 @@ def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): compounds exist, then they will be used. If False and an end group compound is None, then the head or tail port will be exposed in the polymer. - make_periodic : bool, default False - If True and an ``end_groups`` compound is None, then the head and tail - will be forced into an overlap with a periodicity along the z-axis. - If end group compounds exist, then there will be a warning. However - ``add_hydrogens`` will simply be overwritten. - If False, ``end_groups`` compound is None, and ``add_hydrogens`` is + periodic_axis : str, default None + If not ``None`` and an ``end_groups`` compound is None, then the head + and tail will be forced into an overlap with a periodicity along the + ``axis`` (default="z") specified. See + :meth:`mbuild.lib.recipes.polymer.Polymer.create_periodic_bond` for + more details. If end group compounds exist, then there will be a + warning. However ``add_hydrogens`` will simply be overwritten. + If ``None``, ``end_groups`` compound is None, and ``add_hydrogens`` is False then the head or tail port will be exposed in the polymer. """ if n < 1: @@ -194,38 +201,8 @@ def build(self, n, sequence="A", add_hydrogens=True, make_periodic=False): head_tail = [head_port, tail_port] - if make_periodic: - if np.any([x is not None for x in self._end_groups]): - raise ValueError( - "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." - ) - z_axis_transform( - compound=self, - point_on_z_axis=head_port.pos, - point_on_zx_plane=tail_port.pos, - ) - bond_length = 2 * np.max( - [ - np.sqrt(np.sum(np.square(atm1.pos - atm2.pos))) - for (atm1, atm2) in self.root.bond_graph.edges( - self["monomer[0]"] - ) - ] - ) - self.periodic = [ - x > bond_length for x in head_port.pos - tail_port.pos - ] - force_overlap(self, head_port, tail_port) - if not np.all( - np.isclose(self.periodic, [0, 0, 1], np.finfo(float).eps) - ): - raise ValueError( - "Periodic polymer could not be aligned:\n Head Coord {}\n Tail Coord {}\n Difference {}".format( - head_port.pos, - tail_port.pos, - tail_port.pos - head_port.pos, - ) - ) + if periodic_axis is not None: + self.create_periodic_bond(axis=periodic_axis) else: for i, compound in enumerate(self._end_groups): if compound is not None: @@ -411,6 +388,54 @@ def add_end_groups( else: raise ValueError("Label must be 'head' or 'tail'") + def create_periodic_bond(self, axis="z"): + """Align the polymer along a specific axis and connect its end groups. + + If ``end_groups`` is ``None``, then the head and tail will be forced + into an overlap with a periodicity along the ``axis`` specified. If + end group compounds exist, then there will be a warning. However + ``add_hydrogens`` will simply be overwritten. If ``None``, + ``end_groups`` compound is None, and ``add_hydrogens`` is + False then the head or tail port will be exposed in the polymer. + + Parameters + ---------- + axis : str, default="z" + Axis along which to orient the polymer taken as the line connected the + free ports of the end group. May be "x", "y", or "z". + """ + if np.any([x is not None for x in self._end_groups]): + raise ValueError( + "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." + ) + head_port = self.monomers[0]["up"] # First monomer + tail_port = self.monomers[-1]["down"] # Last monomer + if axis == "x": + x_axis_transform( + compound=self, + point_on_x_axis=head_port.pos, + point_on_xy_plane=tail_port.pos, + ) + self.periodicity = (True, self.periodicity[1], self.periodicity[2]) + elif axis == "y": + y_axis_transform( + compound=self, + point_on_y_axis=head_port.pos, + point_on_yx_plane=tail_port.pos, + ) + self.periodicity = (self.periodicity[0], True, self.periodicity[2]) + elif axis == "z": + z_axis_transform( + compound=self, + point_on_z_axis=head_port.pos, + point_on_zx_plane=tail_port.pos, + ) + self.periodicity = (self.periodicity[0], self.periodicity[1], True) + else: + raise ValueError("axis must be either: 'x', 'y', or 'z'") + + force_overlap(self, head_port, tail_port) + def _add_port(compound, label, idx, separation, orientation=None, replace=True): """Add the ports to the compound at the specified particle index. From 433b7458527ed6cb3fae871ba8f17f3cc5d60a4b Mon Sep 17 00:00:00 2001 From: jac16 Date: Fri, 15 Mar 2024 10:58:19 -0400 Subject: [PATCH 07/11] Update to allow choice in periodic axis with a Polymer method --- mbuild/lib/recipes/polymer.py | 29 +++++++++++++++-------------- mbuild/tests/test_monolayer.py | 6 +++--- 2 files changed, 18 insertions(+), 17 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 3918792a9..ca3868b15 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -190,6 +190,7 @@ def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): # Add the end groups head = self["monomer[0]"] # First monomer tail = self["monomer[{}]".format(n_monomers - 1)] # Last monomer + if not head["up"].used: head_port = head["up"] else: @@ -202,7 +203,7 @@ def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): head_tail = [head_port, tail_port] if periodic_axis is not None: - self.create_periodic_bond(axis=periodic_axis) + self.create_periodic_bond(head_port, tail_port, axis=periodic_axis) else: for i, compound in enumerate(self._end_groups): if compound is not None: @@ -388,7 +389,7 @@ def add_end_groups( else: raise ValueError("Label must be 'head' or 'tail'") - def create_periodic_bond(self, axis="z"): + def create_periodic_bond(self, port1, port2, axis="z"): """Align the polymer along a specific axis and connect its end groups. If ``end_groups`` is ``None``, then the head and tail will be forced @@ -400,6 +401,10 @@ def create_periodic_bond(self, axis="z"): Parameters ---------- + port1 : np.array + First port to align the polymer + port2 : np.array + Second port to align the polymer axis : str, default="z" Axis along which to orient the polymer taken as the line connected the free ports of the end group. May be "x", "y", or "z". @@ -408,33 +413,29 @@ def create_periodic_bond(self, axis="z"): raise ValueError( "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." ) - head_port = self.monomers[0]["up"] # First monomer - tail_port = self.monomers[-1]["down"] # Last monomer + if axis == "x": x_axis_transform( compound=self, - point_on_x_axis=head_port.pos, - point_on_xy_plane=tail_port.pos, + new_origin=port1.pos, + point_on_x_axis=port2.pos, ) - self.periodicity = (True, self.periodicity[1], self.periodicity[2]) elif axis == "y": y_axis_transform( compound=self, - point_on_y_axis=head_port.pos, - point_on_yx_plane=tail_port.pos, + new_origin=port1.pos, + point_on_y_axis=port2.pos, ) - self.periodicity = (self.periodicity[0], True, self.periodicity[2]) elif axis == "z": z_axis_transform( compound=self, - point_on_z_axis=head_port.pos, - point_on_zx_plane=tail_port.pos, + new_origin=port1.pos, + point_on_z_axis=port2.pos, ) - self.periodicity = (self.periodicity[0], self.periodicity[1], True) else: raise ValueError("axis must be either: 'x', 'y', or 'z'") - force_overlap(self, head_port, tail_port) + force_overlap(self, port1, port2) def _add_port(compound, label, idx, separation, orientation=None, replace=True): diff --git a/mbuild/tests/test_monolayer.py b/mbuild/tests/test_monolayer.py index 214e37f85..ae1fd5585 100644 --- a/mbuild/tests/test_monolayer.py +++ b/mbuild/tests/test_monolayer.py @@ -50,7 +50,7 @@ def test_pattern_kwargs(self, ch2): def test_periodic_pattern(self, ch2): # Make Periodic without Hydrogen Conflict chain = mb.recipes.Polymer(monomers=[ch2]) - chain.build(n=10, add_hydrogens=False, make_periodic=True) + chain.build(n=10, add_hydrogens=False, periodic_axis="x") assert not chain.all_ports() bonded_atoms = [ @@ -61,7 +61,7 @@ def test_periodic_pattern(self, ch2): # Make Periodic with Hydrogen Conflict chain2 = mb.recipes.Polymer(monomers=[ch2]) - chain2.build(n=10, add_hydrogens=True, make_periodic=True) + chain2.build(n=10, add_hydrogens=True, periodic_axis="y") assert not chain2.all_ports() bonded_atoms = [ @@ -72,7 +72,7 @@ def test_periodic_pattern(self, ch2): with pytest.raises(ValueError): chain3 = mb.recipes.Polymer(monomers=[ch2], end_groups=ch2) - chain3.build(n=10, add_hydrogens=True, make_periodic=True) + chain3.build(n=10, add_hydrogens=True, periodic_axis="z") def test_mixed_monolayer(self, ch2): n = 8 From 6cfea1e3b7a8776d7b3b580deebda079dc062592 Mon Sep 17 00:00:00 2001 From: jac16 Date: Mon, 25 Mar 2024 14:55:58 -0400 Subject: [PATCH 08/11] Expose polymer ports and update tests --- mbuild/lib/recipes/polymer.py | 125 +++++++++++++++------------------ mbuild/tests/test_monolayer.py | 22 +++--- 2 files changed, 66 insertions(+), 81 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index ca3868b15..f1686bd38 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -100,6 +100,8 @@ def __init__(self, monomers=None, end_groups=None): ) self._port_labels = ["up", "down"] self._headtail = [None, None] + self.head_port = None + self.tail_port = None @property def monomers(self): @@ -117,12 +119,16 @@ def end_groups(self): """ return self._end_groups - def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): + def build(self, n, sequence="A", add_hydrogens=True): """Connect one or more components in a specified sequence. Uses the compounds that are stored in Polymer.monomers and Polymer.end_groups. + If no capping method is used, i.e., if ``add_hydrogens == False`` + and ``Polymer.end_groups is None``, then the ports are exposed as, + ``Polymer.head_port`` and ``Polymer.tail_port``. + Parameters ---------- n : int @@ -153,7 +159,6 @@ def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): """ if n < 1: raise ValueError("n must be 1 or more") - n_monomers = n * len(sequence) for monomer in self._monomers: for label in self._port_labels: @@ -187,49 +192,39 @@ def build(self, n, sequence="A", add_hydrogens=True, periodic_axis=None): if n_added == n * len(sequence) - 1: break - # Add the end groups - head = self["monomer[0]"] # First monomer - tail = self["monomer[{}]".format(n_monomers - 1)] # Last monomer - - if not head["up"].used: - head_port = head["up"] - else: - head_port = None - if not tail["down"].used: - tail_port = tail["down"] - else: - tail_port = None - - head_tail = [head_port, tail_port] + self.head_port = first_part["up"] if not first_part["up"].used else None + self.tail_port = ( + last_part["down"] if not last_part["down"].used else None + ) - if periodic_axis is not None: - self.create_periodic_bond(head_port, tail_port, axis=periodic_axis) - else: - for i, compound in enumerate(self._end_groups): - if compound is not None: - if self._headtail[i] is not None: - head_tail[i].update_separation(self._headtail[i]) - self.add(compound) - force_overlap(compound, compound.labels["up"], head_tail[i]) + head_tail = [self.head_port, self.tail_port] + for i, compound in enumerate(self._end_groups): + if compound is not None: + if self._headtail[i] is not None: + head_tail[i].update_separation(self._headtail[i]) + self.add(compound) + force_overlap(compound, compound.labels["up"], head_tail[i]) + head_tail[i] = None + else: + if add_hydrogens: + hydrogen = H() + # Defaut to 1/2 H-C bond len + head_tail[i].update_separation(0.0547) + hydrogen["up"].update_separation(0.0547) + self.add(hydrogen) + force_overlap(hydrogen, hydrogen["up"], head_tail[i]) + head_tail[i] = None else: - if add_hydrogens: - hydrogen = H() - # Defaut to 1/2 H-C bond len - head_tail[i].update_separation(0.0547) - hydrogen["up"].update_separation(0.0547) - self.add(hydrogen) - force_overlap(hydrogen, hydrogen["up"], head_tail[i]) - else: - # if None, hoist port to polymer level - self.add( - head_tail[i], - self._port_labels[i], - containment=False, - ) - - for port in self.all_ports(): - if port not in self.available_ports(): - self.remove(port) + # if None, hoist port to polymer level + self.add( + head_tail[i], + self._port_labels[i], + containment=False, + ) + if head_tail[i] is None and i == 0: + self.head_port = None + elif head_tail[i] is None and i == 1: + self.tail_port = None def add_monomer( self, @@ -258,7 +253,7 @@ def add_monomer( | | H(2) H(5) - If replace=True, then this fucntion removes the hydrogen atoms that are + If replace=True, then this function removes the hydrogen atoms that are occupying where the C-C bond should occur between monomers. It is required that you specify which atoms should be removed which is achieved by the `indices` parameter. @@ -273,7 +268,7 @@ def add_monomer( compound : mbuild.Compound A compound of the individual monomer indices : list of int of length 2 - The particle indicies of compound that represent the polymer + The particle indices of compound that represent the polymer bonding sites. You can specify the indices of particles that will be replaced by the polymer bond, or indices of particles that act as the bonding sites. See the 'replace' parameter notes. @@ -389,53 +384,43 @@ def add_end_groups( else: raise ValueError("Label must be 'head' or 'tail'") - def create_periodic_bond(self, port1, port2, axis="z"): - """Align the polymer along a specific axis and connect its end groups. - - If ``end_groups`` is ``None``, then the head and tail will be forced - into an overlap with a periodicity along the ``axis`` specified. If - end group compounds exist, then there will be a warning. However - ``add_hydrogens`` will simply be overwritten. If ``None``, - ``end_groups`` compound is None, and ``add_hydrogens`` is - False then the head or tail port will be exposed in the polymer. + def create_periodic_bond(self, axis="z"): + """Align and bond the end points of a polymer along an axis. Parameters ---------- - port1 : np.array - First port to align the polymer - port2 : np.array - Second port to align the polymer axis : str, default="z" Axis along which to orient the polymer taken as the line connected the free ports of the end group. May be "x", "y", or "z". """ - if np.any([x is not None for x in self._end_groups]): + if self.head_port is None or self.tail_port is None: raise ValueError( - "Keyword 'make_periodic' cannot be defined if 'end_groups' are provided." + "Polymer head_port and/or tail_port could not be found. Be sure to initialize" + "this object without end_groups and build with add_hydrogens=False." ) - if axis == "x": + if axis.lower() == "x": x_axis_transform( compound=self, - new_origin=port1.pos, - point_on_x_axis=port2.pos, + new_origin=self.head_port.pos, + point_on_x_axis=self.tail_port.pos, ) - elif axis == "y": + elif axis.lower() == "y": y_axis_transform( compound=self, - new_origin=port1.pos, - point_on_y_axis=port2.pos, + new_origin=self.head_port.pos, + point_on_y_axis=self.tail_port.pos, ) - elif axis == "z": + elif axis.lower() == "z": z_axis_transform( compound=self, - new_origin=port1.pos, - point_on_z_axis=port2.pos, + new_origin=self.head_port.pos, + point_on_z_axis=self.tail_port.pos, ) else: raise ValueError("axis must be either: 'x', 'y', or 'z'") - force_overlap(self, port1, port2) + force_overlap(self, self.head_port, self.tail_port) def _add_port(compound, label, idx, separation, orientation=None, replace=True): diff --git a/mbuild/tests/test_monolayer.py b/mbuild/tests/test_monolayer.py index ae1fd5585..3bc84f14f 100644 --- a/mbuild/tests/test_monolayer.py +++ b/mbuild/tests/test_monolayer.py @@ -50,7 +50,8 @@ def test_pattern_kwargs(self, ch2): def test_periodic_pattern(self, ch2): # Make Periodic without Hydrogen Conflict chain = mb.recipes.Polymer(monomers=[ch2]) - chain.build(n=10, add_hydrogens=False, periodic_axis="x") + chain.build(n=10, add_hydrogens=False) + chain.create_periodic_bond(axis="x") assert not chain.all_ports() bonded_atoms = [ @@ -61,18 +62,17 @@ def test_periodic_pattern(self, ch2): # Make Periodic with Hydrogen Conflict chain2 = mb.recipes.Polymer(monomers=[ch2]) - chain2.build(n=10, add_hydrogens=True, periodic_axis="y") - assert not chain2.all_ports() - - bonded_atoms = [ - x.name for x in list(chain2["monomer[0]"][0].direct_bonds()) - ] - assert bonded_atoms.count("H") == 2 - assert bonded_atoms.count("C") == 2 + chain2.build(n=10, add_hydrogens=True) + with pytest.raises(ValueError): + chain2.create_periodic_bond(axis="y") + # Make Periodic with End-Group Conflict + chain3 = mb.recipes.Polymer( + monomers=[ch2], end_groups=[mb.clone(ch2), mb.clone(ch2)] + ) + chain3.build(n=10) with pytest.raises(ValueError): - chain3 = mb.recipes.Polymer(monomers=[ch2], end_groups=ch2) - chain3.build(n=10, add_hydrogens=True, periodic_axis="z") + chain3.create_periodic_bond(axis="z") def test_mixed_monolayer(self, ch2): n = 8 From 64dbb4fe9c079a33ed40e3480bc82d81ec63b6e5 Mon Sep 17 00:00:00 2001 From: jac16 Date: Wed, 3 Apr 2024 17:02:13 -0400 Subject: [PATCH 09/11] Fix Error string --- mbuild/lib/recipes/polymer.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index f1686bd38..92973b051 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -395,7 +395,7 @@ def create_periodic_bond(self, axis="z"): """ if self.head_port is None or self.tail_port is None: raise ValueError( - "Polymer head_port and/or tail_port could not be found. Be sure to initialize" + "Polymer head_port and/or tail_port could not be found. Be sure to initialize " "this object without end_groups and build with add_hydrogens=False." ) From 390e141f8393e412238c69c2498a9b45bf63d8ed Mon Sep 17 00:00:00 2001 From: jac16 Date: Wed, 3 Apr 2024 17:27:52 -0400 Subject: [PATCH 10/11] Update tests --- mbuild/tests/test_monolayer.py | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) diff --git a/mbuild/tests/test_monolayer.py b/mbuild/tests/test_monolayer.py index 3bc84f14f..c27845508 100644 --- a/mbuild/tests/test_monolayer.py +++ b/mbuild/tests/test_monolayer.py @@ -49,10 +49,11 @@ def test_pattern_kwargs(self, ch2): def test_periodic_pattern(self, ch2): # Make Periodic without Hydrogen Conflict - chain = mb.recipes.Polymer(monomers=[ch2]) - chain.build(n=10, add_hydrogens=False) - chain.create_periodic_bond(axis="x") - assert not chain.all_ports() + for axis in ["x", "y", "z"]: + chain = mb.recipes.Polymer(monomers=[ch2]) + chain.build(n=10, add_hydrogens=False) + chain.create_periodic_bond(axis=axis) + assert not chain.all_ports() bonded_atoms = [ x.name for x in list(chain["monomer[0]"][0].direct_bonds()) @@ -74,6 +75,12 @@ def test_periodic_pattern(self, ch2): with pytest.raises(ValueError): chain3.create_periodic_bond(axis="z") + # Make Periodic with Unsupported Axis + chain2 = mb.recipes.Polymer(monomers=[ch2]) + chain2.build(n=10, add_hydrogens=True) + with pytest.raises(ValueError): + chain2.create_periodic_bond(axis="a") + def test_mixed_monolayer(self, ch2): n = 8 m = 8 From e5a136434207f85625a787f0926d67ede71f16ec Mon Sep 17 00:00:00 2001 From: jac16 Date: Thu, 4 Apr 2024 09:41:53 -0400 Subject: [PATCH 11/11] Resolve setting 'first_part' --- mbuild/lib/recipes/polymer.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/mbuild/lib/recipes/polymer.py b/mbuild/lib/recipes/polymer.py index 2e1c6b739..19340fffa 100644 --- a/mbuild/lib/recipes/polymer.py +++ b/mbuild/lib/recipes/polymer.py @@ -178,7 +178,9 @@ def build(self, n, sequence="A", add_hydrogens=True): for n_added, seq_item in enumerate(it.cycle(sequence)): this_part = clone(seq_map[seq_item]) self.add(this_part, "monomer[$]") - if last_part is not None: + if last_part is None: + first_part = this_part + else: # Transform this part, such that its bottom port is rotated # and translated to the last parts top port. force_overlap( @@ -219,7 +221,7 @@ def build(self, n, sequence="A", add_hydrogens=True): self._port_labels[i], containment=False, ) - + if head_tail[i] is None and i == 0: self.head_port = None elif head_tail[i] is None and i == 1: