Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Streamline angles(), bonds(), dihedrals() #3203

Merged
merged 2 commits into from
May 2, 2021
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
101 changes: 16 additions & 85 deletions package/MDAnalysis/core/topologyobjects.py
Original file line number Diff line number Diff line change
Expand Up @@ -899,21 +899,15 @@ def values(self, **kwargs):
elif self.btype == 'improper':
return self.dihedrals(**kwargs)

def _bondsSlow(self, pbc=False): # pragma: no cover
"""Slow version of bond (numpy implementation)"""
if not self.btype == 'bond':
return TypeError("TopologyGroup is not of type 'bond'")
else:
bond_dist = self._ags[0].positions - self._ags[1].positions
if pbc:
box = self._ags[0].dimensions
# orthogonal and divide by zero check
if (box[6:9] == 90.).all() and not (box[0:3] == 0).any():
bond_dist -= np.rint(bond_dist / box[0:3]) * box[0:3]
else:
raise ValueError("Only orthogonal boxes supported")

return np.array([mdamath.norm(a) for a in bond_dist])
def _calc_connection_values(self, func, *btypes, result=None, pbc=False):
if not any(self.btype == btype for btype in btypes):
strbtype = "' or '".join(btypes)
raise TypeError(f"TopologyGroup is not of type '{strbtype}'")
if not result:
result = np.zeros(len(self), np.float64)
box = None if not pbc else self._ags[0].dimensions
positions = [ag.positions for ag in self._ags]
return func(*positions, box=box, result=result)

def bonds(self, pbc=False, result=None):
"""Calculates the distance between all bonds in this TopologyGroup
Expand All @@ -928,30 +922,8 @@ def bonds(self, pbc=False, result=None):

Uses cython implementation
"""
if not self.btype == 'bond':
raise TypeError("TopologyGroup is not of type 'bond'")
if not result:
result = np.zeros(len(self), np.float64)
if pbc:
return distances.calc_bonds(self._ags[0].positions,
self._ags[1].positions,
box=self._ags[0].dimensions,
result=result)
else:
return distances.calc_bonds(self._ags[0].positions,
self._ags[1].positions,
result=result)

def _anglesSlow(self): # pragma: no cover
"""Slow version of angle (numpy implementation)"""
if not self.btype == 'angle':
raise TypeError("TopologyGroup is not of type 'angle'")

vec1 = self._ags[0].positions - self._ags[1].positions
vec2 = self._ags[2].positions - self._ags[1].positions

angles = np.array([mdamath.angle(a, b) for a, b in zip(vec1, vec2)])
return angles
return self._calc_connection_values(distances.calc_bonds, "bond",
pbc=pbc, result=result)

def angles(self, result=None, pbc=False):
"""Calculates the angle in radians formed between a bond
Expand All @@ -975,34 +947,8 @@ def angles(self, result=None, pbc=False):
Added *pbc* option (default ``False``)

"""
if not self.btype == 'angle':
raise TypeError("TopologyGroup is not of type 'angle'")
if not result:
result = np.zeros(len(self), np.float64)
if pbc:
return distances.calc_angles(self._ags[0].positions,
self._ags[1].positions,
self._ags[2].positions,
box=self._ags[0].dimensions,
result=result)
else:
return distances.calc_angles(self._ags[0].positions,
self._ags[1].positions,
self._ags[2].positions,
result=result)

def _dihedralsSlow(self): # pragma: no cover
"""Slow version of dihedral (numpy implementation)"""
if self.btype not in ['dihedral', 'improper']:
raise TypeError("TopologyGroup is not of type 'dihedral' or "
"'improper'")

ab = self._ags[0].positions - self._ags[1].positions
bc = self._ags[1].positions - self._ags[2].positions
cd = self._ags[2].positions - self._ags[3].positions

return np.array([mdamath.dihedral(a, b, c)
for a, b, c in zip(ab, bc, cd)])
return self._calc_connection_values(distances.calc_angles, "angle",
pbc=pbc, result=result)

def dihedrals(self, result=None, pbc=False):
"""Calculate the dihedral angle in radians for this topology
Expand All @@ -1028,21 +974,6 @@ def dihedrals(self, result=None, pbc=False):
.. versionchanged:: 0.9.0
Added *pbc* option (default ``False``)
"""
if self.btype not in ['dihedral', 'improper']:
raise TypeError("TopologyGroup is not of type 'dihedral' or "
"'improper'")
if not result:
result = np.zeros(len(self), np.float64)
if pbc:
return distances.calc_dihedrals(self._ags[0].positions,
self._ags[1].positions,
self._ags[2].positions,
self._ags[3].positions,
box=self._ags[0].dimensions,
result=result)
else:
return distances.calc_dihedrals(self._ags[0].positions,
self._ags[1].positions,
self._ags[2].positions,
self._ags[3].positions,
result=result)
return self._calc_connection_values(distances.calc_dihedrals,
"dihedral", "improper",
pbc=pbc, result=result)