From 1494574183cabe1ffa7a399d65851a01b295cd5c Mon Sep 17 00:00:00 2001 From: avcopan Date: Thu, 15 Aug 2024 20:28:52 -0500 Subject: [PATCH 1/2] New: Track atoms neighboring the torsion axis Working towards being able to replace `None`s with actual neighbors on conversion to geometry. --- automol/data/rotor.py | 7 +- automol/data/tors.py | 283 ++++++++++++++++++++++------------ automol/geom/_conv.py | 6 +- automol/vmat.py | 15 ++ automol/zmat/__init__.py | 2 + automol/zmat/base/__init__.py | 2 + 6 files changed, 213 insertions(+), 102 deletions(-) diff --git a/automol/data/rotor.py b/automol/data/rotor.py index 0768c762..e005adb9 100644 --- a/automol/data/rotor.py +++ b/automol/data/rotor.py @@ -2,6 +2,7 @@ Rotors contain of one or more Torsion data structures """ + import dataclasses import itertools from typing import Dict, List, Optional, Union @@ -210,7 +211,7 @@ def torsion_symmetries(rotor: Rotor) -> List[int]: return tuple(map(tors.symmetry, torsions(rotor))) -def torsion_grids(rotor: Rotor, increment: float=30 * phycon.DEG2RAD) -> List[Grid]: +def torsion_grids(rotor: Rotor, increment: float = 30 * phycon.DEG2RAD) -> List[Grid]: """Get the coordinate grids for the torsions in a rotor :param rotor: A rotor @@ -353,6 +354,7 @@ def rotors_from_zmatrix( """ gra = zmat.graph(zma, stereo=True, dummy=True) if gra is None else gra lin_keys = graph.linear_atom_keys(gra, dummy=True) + znkeys_dct = zmat.neighbor_keys(zma) rot_bkeys = graph.rotational_bond_keys(gra, lin_keys=lin_keys) @@ -369,6 +371,7 @@ def rotors_from_zmatrix( coo=tor_keys, grps=graph.rotational_groups(gra, *tor_axis), symm=graph.rotational_symmetry_number(gra, *tor_axis, lin_keys=lin_keys), + ngrps=list(map(znkeys_dct.get, tor_axis)), ) tor_lst.append(tor) @@ -566,7 +569,7 @@ def rotors_torsion_symmetries( def rotors_torsion_grids( - rotors: List[Rotor], flat: bool = False, increment: float=30 * phycon.DEG2RAD + rotors: List[Rotor], flat: bool = False, increment: float = 30 * phycon.DEG2RAD ) -> Union[List[Grid], List[List[Grid]]]: """Get the torsion coordinate grids from a list of rotors diff --git a/automol/data/tors.py b/automol/data/tors.py index 1983bb66..19d8809e 100644 --- a/automol/data/tors.py +++ b/automol/data/tors.py @@ -4,7 +4,7 @@ """ import dataclasses -from typing import Any, List, Optional, Tuple +from typing import Any import numpy import yaml @@ -14,47 +14,43 @@ from .. import zmat from ..util import ZmatConv, zmat_conv -Axis = Tuple[int, int] -DihCoord = Tuple[int, int, int, int] -Groups = Tuple[List[int], List[int]] -Grid = List[float] +Axis = tuple[int, int] +DihCoord = tuple[int, int, int, int] +Group = list[int] +Groups = tuple[Group, Group] +Grid = list[float] @dataclasses.dataclass class Torsion: - """Encodes information for a single torsion, which is one component of a rotor + """Encodes information for a single torsion, which is one component of a rotor. :param name: The z-matrix coordinate name - :type name: str :param coordinate: The z-matrix keys defining the torsion coordinate - :type coordinate: DihKey :param groups: The sets of atoms keys defining the rotational groups - :type groups: Tuple[List[int], List[int]] - :type symmetry: The rotational symmetry number of the torsion - :type symmetry: int + :param symmetry: The rotational symmetry number of the torsion + :param neighbor_groups: The subsets of rotational group atoms that neighbor the axis """ name: str coordinate: DihCoord groups: Groups symmetry: int + neighbor_groups: Groups | None = None # Torsion functions # # Constructors -def from_data(name_: str, coo: DihCoord, grps: Groups, symm: int) -> Torsion: - """Construct a torsion from data +def from_data( + name_: str, coo: DihCoord, grps: Groups, symm: int, ngrps: Groups | None = None +) -> Torsion: + """Construct a torsion from data. :param name_: The z-matrix coordinate name - :type name_: str :param coo: The z-matrix keys defining the torsion coordinate - :type coo: DihKey :param groups: The sets of atoms keys defining the rotational groups - :type groups: Tuple[List[int], List[int]] :type symm: The rotational symmetry number of the torsion - :type symm: int :return: The torsion data structure - :rtype: Torsion """ assert len(coo) == 4, f"Invalid torsion coordinate: {coo}" assert len(grps) == 2, f"Invalid torsion groups: {grps}" @@ -63,34 +59,29 @@ def from_data(name_: str, coo: DihCoord, grps: Groups, symm: int) -> Torsion: coordinate=tuple(coo), groups=tuple(map(tuple, grps)), symmetry=int(symm), + neighbor_groups=None if ngrps is None else tuple(map(tuple, ngrps)), ) # # Getters def name(tor: Torsion) -> str: - """Get the coordinate name of a torsion + """Get the coordinate name of a torsion. :param tor: A torsion - :type tor: Torsion :return: The coordinate name - :rtype: str """ return tor.name def coordinate( - tor: Torsion, key_typ: str = "zmat", zc_: Optional[ZmatConv] = None + tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None ) -> DihCoord: - """Get the torsion coordinate keys + """Get the torsion coordinate keys. :param tor: A torsion - :type tor: Torsion :param key_typ: The type of keys to return, "zmat" (default) or "geom" - :type key_typ: str, optional :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None - :type zc_: Optional[ZmatConv], optional :return: The torsion rotational keys - :rtype: str """ coo = tor.coordinate if key_typ == "geom": @@ -98,19 +89,13 @@ def coordinate( return coo -def groups( - tor: Torsion, key_typ: str = "zmat", zc_: Optional[ZmatConv] = None -) -> Groups: - """Get the rotational groups of a torsion +def groups(tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None) -> Groups: + """Get the rotational groups of a torsion. :param tor: A torsion - :type tor: Torsion :param key_typ: The type of keys to return, "zmat" (default) or "geom" - :type key_typ: str, optional :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None - :type zc_: Optional[ZmatConv], optional :return: The torsion rotational groups - :rtype: str """ grps = tor.groups if key_typ == "geom": @@ -119,100 +104,140 @@ def groups( def symmetry(tor: Torsion) -> int: - """Get the rotational symmetry of a torsion + """Get the rotational symmetry of a torsion. :param tor: A torsion - :type tor: Torsion :return: The rotational symmetry - :rtype: int """ return tor.symmetry +def neighbor_groups( + tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None +) -> Groups: + """Get the neighbor groups of a torsion. + + :param tor: A torsion + :param key_typ: The type of keys to return, "zmat" (default) or "geom" + :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None + :return: The neighbor groups + """ + ngrps = tor.neighbor_groups + if key_typ == "geom": + ngrps = zmat_conv.relabel_zmatrix_key_sequence(zc_, ngrps) + return ngrps + + # setters def set_name(tor: Torsion, name_: str) -> Torsion: - """Set the coordinate name of a torsion + """Set the coordinate name of a torsion. :param tor: A torsion :param name_: The coordinate name :return: The torsion """ return from_data( - name_=name_, coo=coordinate(tor), grps=groups(tor), symm=symmetry(tor) + name_=name_, + coo=coordinate(tor), + grps=groups(tor), + symm=symmetry(tor), + ngrps=neighbor_groups(tor), ) def set_coordinate(tor: Torsion, coo: DihCoord) -> Torsion: - """get the torsion coordinate keys + """Set the torsion coordinate keys. :param tor: A torsion :param coo: The torsion coordinate keys :return: The torsion """ - return from_data(name_=name(tor), coo=coo, grps=groups(tor), symm=symmetry(tor)) + return from_data( + name_=name(tor), + coo=coo, + grps=groups(tor), + symm=symmetry(tor), + ngrps=neighbor_groups(tor), + ) def set_groups(tor: Torsion, grps: Groups) -> Torsion: - """get the rotational groups of a torsion + """Set the rotational groups of a torsion. :param tor: A torsion :param grps: The torsion rotational groups :return: The torsion """ return from_data( - name_=name(tor), coo=coordinate(tor), grps=grps, symm=symmetry(tor) + name_=name(tor), + coo=coordinate(tor), + grps=grps, + symm=symmetry(tor), + ngrps=neighbor_groups(tor), ) def set_symmetry(tor: Torsion, symm: int) -> Torsion: - """get the rotational symmetry of a torsion + """Set the rotational symmetry of a torsion. :param tor: A torsion :param symm: The rotational symmetry :return: The torsion """ - return from_data(name_=name(tor), coo=coordinate(tor), grps=groups(tor), symm=symm) + return from_data( + name_=name(tor), + coo=coordinate(tor), + grps=groups(tor), + symm=symm, + ngrps=neighbor_groups(tor), + ) + + +def set_neighbor_groups(tor: Torsion, ngrps: Groups) -> Torsion: + """Set the neighbor groups of a torsion. + + :param tor: A torsion + :param ngrps: The neighbor groups + :return: The torsion + """ + return from_data( + name_=name(tor), + coo=coordinate(tor), + grps=groups(tor), + symm=symmetry(tor), + ngrps=ngrps, + ) # properties -def axis(tor: Torsion, key_typ: str = "zmat", zc_: Optional[ZmatConv] = None) -> Axis: - """Get the rotational axis of a torsion +def axis(tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None) -> Axis: + """Get the rotational axis of a torsion. :param tor: A torsion - :type tor: Torsion :param key_typ: The type of keys to return, "zmat" (default) or "geom" - :type key_typ: str, optional :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None - :type zc_: Optional[ZmatConv], optional :return: The torsion rotational axis - :rtype: str """ ks_ = coordinate(tor, key_typ=key_typ, zc_=zc_) return ks_[1:3] def span(tor: Torsion) -> float: - """Get the angular span of a torsion, based on the symmetry number + """Get the angular span of a torsion, based on the symmetry number. :param tor: A torsion - :type tor: Torsion :return: The angular span of the torsion, 2 pi / symmetry - :rtype: float """ return 2 * numpy.pi / symmetry(tor) def grid(tor: Torsion, zma, increment: float = 30 * phycon.DEG2RAD) -> Grid: - """Get the coordinate grid for a torsion + """Get the coordinate grid for a torsion. :param tor: A torsion - :type tor: Torsion :param zma: The z-matrix associated with this torsion - :type zma: automol zmat data structure :param increment: The grid increment, in radians - :type increment: float :return: The coordinate grid - :rtype: Grid """ # [0, 30, 60, 90, ...] << in degrees vals = numpy.arange(0, span(tor), increment) @@ -224,27 +249,25 @@ def grid(tor: Torsion, zma, increment: float = 30 * phycon.DEG2RAD) -> Grid: # # Transformations def with_geometry_indices(tor: Torsion, zc_: ZmatConv) -> Torsion: """Given a z-matrix torsion and a z-matrix conversion, return a torsion with - geometry indices + geometry indices. That is, the axis and group indices will skip dummy atoms :param tor: A torsion data structure - :type tor: Torsion :param zc_: A z-matrix conversion data structure - :type zc_: ZmatConv :return: A torsion data structure using geometry indices - :rtype: Torsion """ return from_data( name_=name(tor), coo=coordinate(tor, key_typ="geom", zc_=zc_), grps=groups(tor, key_typ="geom", zc_=zc_), symm=symmetry(tor), + ngrps=neighbor_groups(tor, key_typ="geom", zc_=zc_), ) def update_against_zmatrix(tor: Torsion, zma: Any) -> Torsion: - """Update a torsion object from a z-matrix + """Update a torsion object from a z-matrix. Temporarily needed to make sure the torsion object contains sufficient information @@ -255,49 +278,45 @@ def update_against_zmatrix(tor: Torsion, zma: Any) -> Torsion: coo = list(reversed(zmat.coordinate(zma, name(tor)))) # In case the coordinate was re-ordered, make sure the groups match grps = tuple(sorted(groups(tor), key=lambda g: coo[-1] in g)) - return from_data(name_=name(tor), coo=coo, grps=grps, symm=symmetry(tor)) + # Re-evaluate the neighbor groups, in case they are missing + nkeys_dct = zmat.neighbor_keys(zma) + ngrps = tuple(sorted(nkeys_dct[k]) for k in coo[1:3]) + return from_data( + name_=name(tor), coo=coo, grps=grps, symm=symmetry(tor), ngrps=ngrps + ) -# Torsion List functions -def torsions_string(tor_lst: List[Torsion], one_indexed: bool = True) -> str: - """Write a list of torsions to a string +# Torsion list functions +def torsions_string(tor_lst: list[Torsion], one_indexed: bool = True) -> str: + """Write a list of torsions to a string. :param tor_lst: A list of torsions - :type tor_lst: List[Torsion] :param one_indexed: Is this a one-indexed string? defaults to True - :type one_indexed: bool, optional :returns: A string representations of the torsions, as a flattened list - :rtype: str """ tor_yml_dct = torsions_yaml_data(tor_lst, one_indexed=one_indexed) tor_str = yaml.dump(tor_yml_dct, sort_keys=False) return tor_str -def torsions_from_string(tor_str: str, one_indexed: bool = True) -> List[Torsion]: - """Write a list of torsions to a string +def torsions_from_string(tor_str: str, one_indexed: bool = True) -> list[Torsion]: + """Write a list of torsions to a string. :param tor_str: A string representations of the torsions, as a flattened list - :type tor_str: str :param one_indexed: Is this a one-indexed string? defaults to True - :type one_indexed: bool, optional :returns: A list of torsions - :rtype: List[Torsion] """ tor_yml_dct = yaml.load(tor_str, Loader=yaml.FullLoader) tor_lst = torsions_from_yaml_data(tor_yml_dct, one_indexed=one_indexed) return tor_lst -def torsions_yaml_data(tor_lst: List[Torsion], one_indexed: bool = True) -> dict: - """Write a list of torsions to a yaml-formatted dictionary +def torsions_yaml_data(tor_lst: list[Torsion], one_indexed: bool = True) -> dict: + """Write a list of torsions to a yaml-formatted dictionary. :param tor_lst: A list of torsions - :type tor_lst: List[Torsion] :param one_indexed: Is this a one-indexed string? defaults to True - :type one_indexed: bool, optional :returns: A string representations of the torsions, as a flattened list - :rtype: str """ shift = 1 if one_indexed else 0 @@ -305,40 +324,39 @@ def torsions_yaml_data(tor_lst: List[Torsion], one_indexed: bool = True) -> dict for tor in tor_lst: coo = [k + shift for k in coordinate(tor)] axis1, axis2 = (k + shift for k in axis(tor)) - grps = ([k + shift for k in g] for g in groups(tor)) - grp1, grp2 = ("-".join(map(str, g)) if len(g) > 1 else g[0] for g in grps) coo_str = "-".join("*" if k is None else str(k) for k in coo) - tor_yml_dct[name(tor)] = { + yml_dct = { "axis1": axis1, - "group1": grp1, "axis2": axis2, - "group2": grp2, "symmetry": symmetry(tor), "coordinate": coo_str, } + yml_dct = _write_groups_to_yaml_dict( + yml_dct, groups(tor), shift=shift, key_prefix="group" + ) + yml_dct = _write_groups_to_yaml_dict( + yml_dct, neighbor_groups(tor), shift=shift, key_prefix="neighbor_group" + ) + tor_yml_dct[name(tor)] = yml_dct return tor_yml_dct def torsions_from_yaml_data(tor_yml_dct: dict, one_indexed: bool = True) -> dict: - """Read a list of torsions out of a yaml-formatted torsion dictionary + """Read a list of torsions out of a yaml-formatted torsion dictionary. :param tor_lst: A list of torsions - :type tor_lst: List[Torsion] :param one_indexed: Is this a one-indexed string? defaults to True - :type one_indexed: bool, optional :returns: A string representations of the torsions, as a flattened list - :rtype: str """ shift = -1 if one_indexed else 0 tor_lst = [] - for name_, vals_dct in tor_yml_dct.items(): - ax_ = list(map(vals_dct.__getitem__, ["axis1", "axis2"])) + for name_, yml_dct in tor_yml_dct.items(): + ax_ = list(map(yml_dct.get, ["axis1", "axis2"])) ax_ = [k + shift for k in ax_] - grps = list(map(vals_dct.__getitem__, ["group1", "group2"])) - grps = [[g] if isinstance(g, int) else map(int, g.split("-")) for g in grps] - grps = [[k + shift for k in g] for g in grps] - coo = vals_dct.get("coordinate", None) + grps = _read_groups_from_yaml_dict(yml_dct, shift, "group") + ngrps = _read_groups_from_yaml_dict(yml_dct, shift, "neighbor_group") + coo = yml_dct.get("coordinate", None) if coo is None: coo = [None, *ax_, None] else: @@ -348,9 +366,78 @@ def torsions_from_yaml_data(tor_yml_dct: dict, one_indexed: bool = True) -> dict name_=name_, coo=coo, grps=grps, - symm=vals_dct["symmetry"], + symm=yml_dct["symmetry"], + ngrps=ngrps, ) tor_lst.append(tor) return tuple(tor_lst) + + +# helper functions +def _write_groups_to_yaml_dict( + yml_dct: dict[object, object], + grps: Groups | None, + shift: int, + key_prefix: str, +) -> dict[object, object]: + """Write rotational groups to a YAML dictionary. + + :param yml_dct: YAML dictionary + :param grps: Rotational groups + :param shift: The shift to apply to the keys + :param key_prefix: Prefix for key + :return: Updated YAML dictionary + """ + if grps is None: + return yml_dct + + yml_dct = yml_dct.copy() + grp_vals = tuple(_group_to_yaml_value(g, shift) for g in grps) + yml_dct.update({f"{key_prefix}{i+1}": v for i, v in enumerate(grp_vals)}) + return yml_dct + + +def _read_groups_from_yaml_dict( + yml_dct: dict[object, object], + shift: int, + key_prefix: str, +) -> Groups | None: + """Read rotational groups from a YAML dictionary. + + :param yml_dct: YAML dictionary + :param shift: The shift to apply to the keys + :param key_prefix: Prefix for key + :return: Updated YAML dictionary + """ + grp_vals = [yml_dct.get(f"{key_prefix}{i+1}") for i in range(2)] + any_none = any(v is None for v in grp_vals) + all_none = all(v is None for v in grp_vals) + assert all_none or not any_none, f"In consistent group values: {yml_dct}" + if any_none: + return None + return tuple(_group_from_yaml_value(v, shift) for v in grp_vals) + + +def _group_to_yaml_value(grp: Group, shift: int) -> str | int: + """Get a YAML value to represent a group. + + :param grp: The rotational group + :param shift: The shift to apply to the keys + :return: The YAML value for the group + """ + grp = [k + shift for k in grp] + return "-".join(map(str, grp)) if len(grp) > 1 else grp[0] + + +def _group_from_yaml_value(grp_val: str | int, shift: int) -> Group: + """Get a group from a YAML value. + + :param grp_str: The string representation of the rotational group + :param shift: The shift to apply to the keys + :return: The rotational group + """ + grp = [grp_val] if isinstance(grp_val, int) else map(int, grp_val.split("-")) + grp = [k + shift for k in grp] + return grp diff --git a/automol/geom/_conv.py b/automol/geom/_conv.py index 1f3deea4..a0cf663d 100644 --- a/automol/geom/_conv.py +++ b/automol/geom/_conv.py @@ -624,8 +624,10 @@ def display( # If requested, visualize only a subset of the bonds by removing others if vis_bkeys is not None: - excl_bkeys = graph_base.bond_keys(gra) - set(map(frozenset, vis_bkeys)) - gra = graph_base.remove_bonds(gra, excl_bkeys, stereo=False, check=False) + all_bkeys = graph_base.bond_keys(gra) + vis_bkeys = set(map(frozenset, vis_bkeys)) + gra = graph_base.remove_bonds(gra, all_bkeys - vis_bkeys, stereo=False) + gra = graph_base.add_bonds(gra, vis_bkeys - all_bkeys) view = py3dmol_view(geo, gra=gra, view=view, image_size=image_size, mode=mode) diff --git a/automol/vmat.py b/automol/vmat.py index f83fb876..0f04ea17 100644 --- a/automol/vmat.py +++ b/automol/vmat.py @@ -1,6 +1,7 @@ """V-Matrix: Variable V-Matrix (V-Matrix without coordinate values).""" import itertools +from collections import defaultdict import more_itertools import numpy @@ -524,6 +525,20 @@ def conversion_info(zma: VMatrix) -> ZmatConv: return zmat_conv.from_zmat_data(zcount, src_zkeys_dct) +def neighbor_keys(vma: VMatrix) -> dict[Key, frozenset[Key]]: + """Identify which atoms are explicit neighbors in the V-Matrix. + + :param vma: V-Matrix + :return: A dictionary mapping atoms onto their neighbors + """ + dist_coos = list(distance_coordinates(vma).values()) + nkeys_dct = defaultdict(set) + for key1, key2 in dist_coos: + nkeys_dct[key1].add(key2) + nkeys_dct[key2].add(key1) + return dict_.transform_values(nkeys_dct, frozenset) + + # # V-Matrix-specific functions # # # setters def set_key_matrix(vma: VMatrix, key_mat: KeyMatrix) -> VMatrix: diff --git a/automol/zmat/__init__.py b/automol/zmat/__init__.py index 54d091b6..9db6a2bc 100644 --- a/automol/zmat/__init__.py +++ b/automol/zmat/__init__.py @@ -39,6 +39,7 @@ from ..vmat import dummy_keys from ..vmat import dummy_source_dict from ..vmat import conversion_info +from ..vmat import neighbor_keys # core functions # # constructors from .base._core import from_data @@ -172,6 +173,7 @@ 'dummy_keys', 'dummy_source_dict', 'conversion_info', + 'neighbor_keys', # extra base functions 'samples', 'constraint_dict', diff --git a/automol/zmat/base/__init__.py b/automol/zmat/base/__init__.py index 06b5aca2..9d17f133 100644 --- a/automol/zmat/base/__init__.py +++ b/automol/zmat/base/__init__.py @@ -32,6 +32,7 @@ from ...vmat import dummy_keys from ...vmat import dummy_source_dict from ...vmat import conversion_info +from ...vmat import neighbor_keys # core functions # # constructors from ._core import from_data @@ -136,6 +137,7 @@ 'dummy_keys', 'dummy_source_dict', 'conversion_info', + 'neighbor_keys', # extra functions 'samples', 'constraint_dict', From 0f4c17051859088b9ca363b81c5ae599b868bc8f Mon Sep 17 00:00:00 2001 From: avcopan Date: Thu, 15 Aug 2024 21:30:20 -0500 Subject: [PATCH 2/2] Fix: Replace torsion dummies on zmat->geom conversion Previously, we were putting in `None`, which was causing problems. --- automol/data/tors.py | 22 ++++++++++++++++++++-- automol/tests/test_rotor.py | 14 +++++++++++++- 2 files changed, 33 insertions(+), 3 deletions(-) diff --git a/automol/data/tors.py b/automol/data/tors.py index 19d8809e..a4c24aa6 100644 --- a/automol/data/tors.py +++ b/automol/data/tors.py @@ -74,17 +74,30 @@ def name(tor: Torsion) -> str: def coordinate( - tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None + tor: Torsion, + key_typ: str = "zmat", + zc_: ZmatConv | None = None, + replace_dummy: bool = True, ) -> DihCoord: """Get the torsion coordinate keys. :param tor: A torsion :param key_typ: The type of keys to return, "zmat" (default) or "geom" :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None + :param replace_dummy: Replace dummy atom with other (presumably in-line) neighbor? :return: The torsion rotational keys """ coo = tor.coordinate if key_typ == "geom": + if replace_dummy: + dkeys = zmat_conv.dummy_keys(zc_) + ngrps = neighbor_groups(tor, dummy=False, zc_=zc_) + end_keys0 = [coo[0], coo[-1]] + end_keys = [ + k if k not in dkeys else next(iter(nks)) + for k, nks in zip(end_keys0, ngrps, strict=True) + ] + coo = (end_keys[0], coo[1], coo[2], end_keys[1]) coo = zmat_conv.relabel_zmatrix_key_sequence(zc_, coo, dummy=True) return coo @@ -113,16 +126,21 @@ def symmetry(tor: Torsion) -> int: def neighbor_groups( - tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None + tor: Torsion, key_typ: str = "zmat", zc_: ZmatConv | None = None, dummy: bool = True ) -> Groups: """Get the neighbor groups of a torsion. :param tor: A torsion :param key_typ: The type of keys to return, "zmat" (default) or "geom" :param zc_: Z-matrix conversion info, to avoid re-calculation, defaults to None + :param dummy: Include dummy atom neighbors? :return: The neighbor groups """ ngrps = tor.neighbor_groups + if not dummy: + dkeys = zmat_conv.dummy_keys(zc_) + ngrps = tuple([k for k in g if k not in dkeys] for g in ngrps) + if key_typ == "geom": ngrps = zmat_conv.relabel_zmatrix_key_sequence(zc_, ngrps) return ngrps diff --git a/automol/tests/test_rotor.py b/automol/tests/test_rotor.py index 97f99dc8..e4218f2a 100644 --- a/automol/tests/test_rotor.py +++ b/automol/tests/test_rotor.py @@ -393,10 +393,22 @@ def test__consistency(): assert tor_grps == ref_tor_grps, f"{tor_grps} != {ref_tor_grps}" +def test__dummy_replacement(): + """Check that dummy atoms get replaced upon conversion to geometry.""" + zma = C2H5OH_CH3_ZMA + zgra = reac.ts_graph(C2H5OH_CH3_ZRXN) + rotors = rotor.rotors_from_zmatrix(zma, gra=zgra) + + tor_coos = rotor.rotors_torsion_coordinates(rotors, key_typ="geom") + ref_tor_coos = (((2, 0, 1, 3),), ((0, 1, 3, 6),), ((3, 6, 9, 10),)) + assert tor_coos == ref_tor_coos, f"{tor_coos} != {ref_tor_coos}" + + if __name__ == "__main__": # test__rotor() # test__rotor_with_dummy_atoms() # test__rotor_multidimensional() # test__torsion_list_string() # test__rotor_for_ts() - test__consistency() + # test__consistency() + test__dummy_replacement()