From ee533379daed7f5d6d62edab1c22a2b69d8ab939 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 16:45:12 +0100 Subject: [PATCH 01/11] feat: add checkBlocks for validating block assignments Adds a checker for the blocks used by the rotating translating blocks models, so that they are validated in one place rather than partly and late. It checks that the blocks are list like, that their identifiers are all of one type, that they are one dimensional, that there is one per node, and that they leave more degrees of freedom than a rigid body has. The identifiers are counted rather than used as indices, so they can be of any hashable type, such as strings. Mixed types are rejected because counting them by value and coercing them into an array disagree: 0 and '0' are two blocks by value but one after coercion. The degrees of freedom are counted the way the block Hessian itself is built, where a block holding a single node contributes 3 rather than 6, having nothing to rotate. Two blocks of one node each therefore leave the same 6 degrees of freedom as a single block of everything, so counting the blocks alone would not catch it. Co-Authored-By: Claude Opus 5 --- prody/tests/utilities/test_checkers.py | 82 +++++++++++++++++++++++++- prody/utilities/checkers.py | 69 +++++++++++++++++++++- 2 files changed, 149 insertions(+), 2 deletions(-) diff --git a/prody/tests/utilities/test_checkers.py b/prody/tests/utilities/test_checkers.py index 5726ba258..d8cb43cb2 100644 --- a/prody/tests/utilities/test_checkers.py +++ b/prody/tests/utilities/test_checkers.py @@ -1,13 +1,17 @@ from prody.tests import TestCase +from numpy import arange, array, zeros from numpy.random import random from numpy.testing import assert_equal -from prody.utilities import checkCoords, checkTypes +from prody.utilities import checkBlocks, checkCoords, checkTypes COORDS = random((10, 3))*10 COORDSET = random((2, 10, 3))*10 +NNODES = 4 +BLOCKS = [0, 0, 1, 1] + class TestCheckCoords(TestCase): def testInvalidCoords(self): @@ -29,6 +33,82 @@ def testCoordsetNatoms(self): natoms=20) +class TestCheckBlocks(TestCase): + + def testBlocksList(self): + + self.assertIsNone(checkBlocks(BLOCKS, NNODES)) + + def testBlocksTuple(self): + + self.assertIsNone(checkBlocks(tuple(BLOCKS), NNODES)) + + def testBlocksArray(self): + + self.assertIsNone(checkBlocks(array(BLOCKS), NNODES)) + + def testBlocksNotNumbers(self): + + self.assertIsNone(checkBlocks(['a', 'a', 'b', 'b'], NNODES)) + + def testBlocksOnePerNode(self): + + self.assertIsNone(checkBlocks(arange(NNODES), NNODES)) + + def testInvalidBlocks(self): + + self.assertRaises(TypeError, checkBlocks, None, NNODES) + + def testScalarBlocks(self): + + self.assertRaises(TypeError, checkBlocks, 0, NNODES) + + def testStringBlocks(self): + + self.assertRaises(TypeError, checkBlocks, 'aabb', NNODES) + + def testNestedBlocks(self): + + self.assertRaises(ValueError, checkBlocks, [[0, 1], [1, 0]], NNODES) + + def testBlocks2D(self): + + self.assertRaises(ValueError, checkBlocks, zeros((NNODES, 2), int), + NNODES) + + def testTooFewBlocks(self): + + self.assertRaises(ValueError, checkBlocks, BLOCKS[:-1], NNODES) + + def testTooManyBlocks(self): + + self.assertRaises(ValueError, checkBlocks, BLOCKS + [1], NNODES) + + def testOneBlock(self): + + self.assertRaises(ValueError, checkBlocks, zeros(NNODES, int), NNODES) + + def testTwoBlocksOfOneNode(self): + # two blocks of one node have 3 degrees of freedom each, which is no + # more than the 6 of a rigid body + + self.assertRaises(ValueError, checkBlocks, [0, 1], 2) + + def testMixedTypeBlocks(self): + + self.assertRaises(TypeError, checkBlocks, [0, 0, '1', '1'], NNODES) + + def testRaggedBlocks(self): + + self.assertRaises(ValueError, checkBlocks, [[0, 1], [2]], 2) + + def testNodesLabel(self): + + self.assertRaisesRegex(ValueError, 'coarse grained nodes', + checkBlocks, BLOCKS[:-1], NNODES, + 'coarse grained nodes') + + class testCheckTypes(TestCase): def testCorrectMonotypeOneArg(self): diff --git a/prody/utilities/checkers.py b/prody/utilities/checkers.py index 16cea72c3..23d13d500 100644 --- a/prody/utilities/checkers.py +++ b/prody/utilities/checkers.py @@ -1,8 +1,18 @@ """This module defines functions for type, value, and/or attribute checking.""" +from collections import Counter + from numpy import any, float32, tile -__all__ = ['checkCoords', 'checkWeights', 'checkTypes', 'checkAnisous'] +from .misctools import isListLike + +__all__ = [ + 'checkAnisous', + 'checkBlocks', + 'checkCoords', + 'checkTypes', + 'checkWeights', +] COORDS_NDIM = set([2]) CSETS_NDIMS = set([2, 3]) @@ -102,6 +112,63 @@ def checkWeights(weights, natoms, ncsets=None, dtype=float): return weights +def checkBlocks(blocks, nnodes, label='atoms'): + """Raises an exception unless *blocks* assigns each of *nnodes* nodes to a + block, leaving more degrees of freedom than a rigid body has, otherwise + returns **None**. + + Blocks that leave only the 6 degrees of freedom of a rigid body, such as a + single block of everything or two blocks of one node each, give the block + based elastic network models no internal motion to describe, so they are + rejected here rather than left to fail later on when the block Hessian is + decomposed. Note that a block of a single node has only 3 of those degrees + of freedom, having nothing to rotate. + + :arg blocks: a list or array of block identifiers, one per node. These are + counted rather than used as indices, so they can be of any hashable + type, such as strings, as long as they are all of the same one. + :type blocks: list, :class:`numpy.ndarray` + + :arg nnodes: the number of nodes *blocks* has to assign + :type nnodes: int + + :arg label: what the nodes are called in the exceptions raised, + default is ``'atoms'`` + :type label: str + """ + + if not isListLike(blocks): + raise TypeError('blocks must be list like, not %s' + % type(blocks).__name__) + + # identifiers of mixed types would be counted separately here but could + # still be coerced into agreeing elsewhere, so they are not allowed + types = set(type(block) for block in blocks) + if len(types) > 1: + raise TypeError('blocks must all be of the same type, not %s' + % ', '.join(sorted(type_.__name__ for type_ in types))) + + # every identifier is of the same type by now, so one stands for them all + if len(blocks) and isListLike(blocks[0]): + raise ValueError('blocks must be one dimensional, giving one identifier ' + 'for each node rather than a list of them') + + if len(blocks) != nnodes: + raise ValueError('len(blocks) must match number of %s, %d, not %d' + % (label, nnodes, len(blocks))) + + # counted as the blocks themselves are, where a block of a single node + # contributes 3 degrees of freedom instead of 6, having nothing to rotate + counter = Counter(blocks) + nsingle = sum(1 for size in counter.values() if size == 1) + ndof = len(counter) * 6 - nsingle * 3 + + if ndof <= 6: + raise ValueError('blocks must leave more than the 6 degrees of freedom ' + 'of a rigid body, not %d, otherwise there is no ' + 'internal motion for them to describe' % ndof) + + def checkTypes(args, **types): """Returns **True** if types of all *args* match those given in *types*. From 36c16db7e0c95cddc37dab2af6abb4ea66201733 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 16:45:12 +0100 Subject: [PATCH 02/11] fix: reject degenerate RTB and imANM blocks before decomposing RTB checked that the number of blocks matched the number of atoms but not that they left anything to sample, so blocks describing a single rigid body were projected and only failed afterwards in calcModes, with an eigenvalue index range that says nothing about the blocks: ValueError: Requested eigenvalue indices are not valid. Valid range is [0, 5] and start <= end, but start=0, end=7 is given Both checks now come from checkBlocks, which imANM inherits, so the blocks are rejected where they are given. Adds the first tests of imANM, covering blocks of the wrong length, blocks that leave a rigid body, identifiers that are not numbers, and identifiers of mixed types. Co-Authored-By: Claude Opus 5 --- prody/dynamics/rtb.py | 5 +- prody/tests/dynamics/test_enms.py | 82 +++++++++++++++++++++++++++++++ 2 files changed, 84 insertions(+), 3 deletions(-) diff --git a/prody/dynamics/rtb.py b/prody/dynamics/rtb.py index 18ba51632..213003977 100644 --- a/prody/dynamics/rtb.py +++ b/prody/dynamics/rtb.py @@ -5,7 +5,7 @@ import numpy as np from prody import LOGGER -from prody.utilities import checkCoords +from prody.utilities import checkBlocks, checkCoords from .anm import ANMBase from .gnm import GNMBase @@ -76,8 +76,7 @@ def buildHessian(self, coords, blocks, cutoff=15., gamma=1., **kwargs): def calcProjection(self, coords, blocks, **kwargs): natoms = self._n_atoms - if natoms != len(blocks): - raise ValueError('len(blocks) must match number of atoms') + checkBlocks(blocks, natoms) LOGGER.timeit('_rtb') from collections import defaultdict diff --git a/prody/tests/dynamics/test_enms.py b/prody/tests/dynamics/test_enms.py index cb9a602b9..efd29c546 100644 --- a/prody/tests/dynamics/test_enms.py +++ b/prody/tests/dynamics/test_enms.py @@ -372,6 +372,88 @@ class TestGNMCalcModes(unittest.TestCase): def setUp(self): pass +class TestRTBBlocks(unittest.TestCase): + """Tests of the blocks given to :meth:`.RTB.buildHessian`, which imANM + inherits. These need no force field, so they are checked here rather than + through ClustRTB and ClustImANM, which cannot reach them without first + fixing a structure.""" + + BLOCK_CLASSES = (RTB, imANM) + + def testBlocksWrongLength(self): + """Test response to blocks that do not cover every atom.""" + for cls in self.BLOCK_CLASSES: + blocks = np.zeros(ATOMS2.numAtoms() - 5, dtype=int) + with self.assertRaisesRegex( + ValueError, 'must match number of atoms', + msg='%s given too few blocks failed to say so' + % cls.__name__): + cls().buildHessian(ATOMS2, blocks) + + def testBlocksAllOneValue(self): + """Test response to blocks that are all the same value. + + One block makes the whole structure a single rigid body, which has no + internal degrees of freedom, so this has to be rejected here rather + than left to fail in the eigendecomposition downstream.""" + for cls in self.BLOCK_CLASSES: + blocks = np.zeros(ATOMS2.numAtoms(), dtype=int) + with self.assertRaisesRegex( + ValueError, 'degrees of freedom of a rigid body', + msg='%s given a single block failed to say so' + % cls.__name__): + cls().buildHessian(ATOMS2, blocks) + + def testBlocksOfOneNode(self): + """Test response to blocks that each hold a single atom. + + Such a block has only 3 degrees of freedom, having nothing to rotate, + so two of them are no better than one block of everything.""" + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'degrees of freedom of a rigid body', + msg='%s given two blocks of one atom failed to say so' + % cls.__name__): + cls().buildHessian(ATOMS2[:2], np.arange(2)) + + def testBlocksNotNumbers(self): + """Test that block identifiers do not have to be numbers""" + for cls in self.BLOCK_CLASSES: + n = ATOMS2.numAtoms() + enm = cls() + enm.buildHessian(ATOMS2, ['a' if i < n // 2 else 'b' + for i in range(n)]) + assert_equal(enm._dof, 12, + '%s with two named blocks failed to give the degrees ' + 'of freedom of two blocks' % cls.__name__) + + def testBlocksMixedTypes(self): + """Test response to block identifiers of more than one type.""" + for cls in self.BLOCK_CLASSES: + n = ATOMS2.numAtoms() + blocks = [0 if i < n // 2 else 'b' for i in range(n)] + with self.assertRaises( + TypeError, + msg='%s given blocks of mixed types failed to raise a ' + 'TypeError' % cls.__name__): + cls().buildHessian(ATOMS2, blocks) + + def testBlocksTwoValues(self): + """Test that two blocks give a projection to sample along""" + for cls in self.BLOCK_CLASSES: + n = ATOMS2.numAtoms() + enm = cls() + enm.buildHessian(ATOMS2, np.arange(n) // (n // 2)) + + projection = enm._getProjection() + assert_equal(projection.shape[0], n * 3, + '%s with two blocks failed to project every degree ' + 'of freedom' % cls.__name__) + self.assertGreater(projection.shape[1], 0, + '%s with two blocks failed to leave any block ' + 'degrees of freedom' % cls.__name__) + + class TestRTB(unittest.TestCase): def testHessian(self): From 37703ec207d5e2a8673373dc5a683c29fa728325 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 16:46:03 +0100 Subject: [PATCH 03/11] fix: check ClustENM inputs and name the fixed structure it writes writePDBFixed now takes the filename to write, and returns the one it wrote, so that a caller no longer has to reconstruct the name from the title. runANMD needs this: it replaces the spaces of the title with underscores while writePDBFixed did not, so a title containing a space made the two disagree and the run failed with FileNotFoundError on a file that had just been written under another name. The file it writes is also closed now rather than left to the garbage collector. Three input checks in run(): - the threshold size mismatch formatted %d against the tuple of thresholds instead of its length, so it raised TypeError from inside its own error message rather than the ValueError it meant to. The maxclust branch above it already used len(). - running without having set atoms gave TypeError: 'NoneType' object is not subscriptable from indexing them for the Kirchhoff matrix. - the blocks of the block based subclasses were only checked once their ANM was built, after a structure had been fixed and a generation minimised. They are now checked against the coarse grained nodes they have to describe before any of that, and the count they are checked against is named in the message, being the nodes rather than every atom of the fixed structure. Adds tests for ClustENM, covering run() arguments, the state of a new instance, titles, the block based subclasses, writePDBFixed, fixing a structure, and a small run. The subclass checks for blocks are made again with blocks set, where the parent's checks would otherwise be masked. Zero generations is covered as the supported way to minimise a structure without sampling, which is why the clustering parameters are not needed in that case. Co-Authored-By: Claude Opus 5 --- prody/dynamics/clustenm.py | 37 +- prody/tests/dynamics/test_clustenm.py | 729 ++++++++++++++++++++++++++ 2 files changed, 758 insertions(+), 8 deletions(-) create mode 100644 prody/tests/dynamics/test_clustenm.py diff --git a/prody/dynamics/clustenm.py b/prody/dynamics/clustenm.py index f5d5c712e..6094fa2b6 100644 --- a/prody/dynamics/clustenm.py +++ b/prody/dynamics/clustenm.py @@ -47,7 +47,7 @@ from prody.measure import calcTransformation, applyTransformation, calcRMSD from prody.ensemble import Ensemble from prody.proteins import writePDB, parsePDB, writePDBStream, parsePDBStream -from prody.utilities import createStringIO, importLA, mad +from prody.utilities import checkBlocks, createStringIO, importLA, mad la = importLA() norm = la.norm @@ -950,19 +950,27 @@ def _getCoordsets(self, indices=None, selected=True): return super(ClustENM, self)._getCoordsets(I, selected) - def writePDBFixed(self): + def writePDBFixed(self, filename=None, replace_spaces=False): 'Write the fixed (initial) structure to a pdb file.' + if filename is None: + filename = self.getTitle()[:-8] + 'fixed.pdb' + if replace_spaces: + filename = filename.replace(" ", "_") + try: from openmm.app import PDBFile except ImportError: raise ImportError('Please install PDBFixer and OpenMM 7.6 in order to use ClustENM.') - PDBFile.writeFile(self._topology, - self._positions, - open(self.getTitle()[:-8] + 'fixed.pdb', 'w'), - keepIds=True) + with open(filename, 'w') as stream: + PDBFile.writeFile(self._topology, + self._positions, + stream, + keepIds=True) + + return filename def writePDB(self, filename=None, single=True, **kwargs): @@ -1186,7 +1194,7 @@ def run(self, cutoff=15., n_modes=3, gamma=1., n_confs=50, rmsd=1.0, if len(self._maxclust) != self._n_gens + 1: raise ValueError( - 'size mismatch: %d generations were set; %d maxclusts were given' % ( + 'size mismatch: %d generations were set; %d maxclusts were given' % ( # noqa: UP031 self._n_gens + 1, len(self._maxclust))) if threshold is None: @@ -1198,7 +1206,9 @@ def run(self, cutoff=15., n_modes=3, gamma=1., n_confs=50, rmsd=1.0, self._threshold = (0,) + (threshold,) * n_gens if len(self._threshold) != self._n_gens + 1: - raise ValueError('size mismatch: %d generations were set; %d thresholds were given' % (self._n_gens + 1, self._threshold)) + raise ValueError( + 'size mismatch: %d generations were set; %d thresholds were given' % ( # noqa: UP031 + self._n_gens + 1, len(self._threshold))) self._sol = solvent if self._nuc is None else 'exp' self._padding = kwargs.pop('padding', 1.0) @@ -1225,6 +1235,17 @@ def run(self, cutoff=15., n_modes=3, gamma=1., n_confs=50, rmsd=1.0, self._cycle = 0 + if self._atoms is None: + raise ValueError('atoms are not set, use `setAtoms`') + + # the block based subclasses check that blocks are set before getting + # here, so this checks them against the nodes they have to describe, + # which fails the run before it minimises anything rather than later + # when its ANM is built + blocks = getattr(self, '_blocks', None) + if blocks is not None: + checkBlocks(blocks, self._n_cg, 'coarse grained nodes') + # check for discontinuity in the structure gnm = GNM() gnm.buildKirchhoff(self._atoms[self._idx_cg], cutoff=self._cutoff) diff --git a/prody/tests/dynamics/test_clustenm.py b/prody/tests/dynamics/test_clustenm.py new file mode 100644 index 000000000..237593946 --- /dev/null +++ b/prody/tests/dynamics/test_clustenm.py @@ -0,0 +1,729 @@ +"""This module contains unit tests for :mod:`~prody.dynamics.clustenm`.""" + +import os +import shutil +import sys +import tempfile +import types +import warnings +import numpy as np +from numpy.testing import * +from prody.utilities import importDec +dec = importDec() + +import prody +from prody import * +from prody import LOGGER +from prody.dynamics.clustenm import ClustENM, ClustRTB, ClustImANM +from prody.tests import unittest +from prody.tests.datafiles import * + +# Import mock tools +try: + from unittest.mock import MagicMock, patch +except ImportError: + from mock import MagicMock, patch + +# Prevent threading hangs on remote servers +os.environ['OMP_NUM_THREADS'] = '1' +os.environ['MKL_NUM_THREADS'] = '1' + +LOGGER.verbosity = 'none' + +# The title a ClustENM gets from setAtoms is the title of the atoms plus +# '_clustenm', and writePDBFixed strips those 8 characters back off to build +# its default filename. Setting the title directly gives the same result +# without needing PDBFixer to fix a real structure. +TITLE = 'my prot' +CLUSTENM_TITLE = TITLE + '_clustenm' +SPACED_NAME = TITLE + '_fixed.pdb' +UNDERSCORED_NAME = TITLE.replace(' ', '_') + '_fixed.pdb' + +# Fixing and minimising a structure needs both of these. They are optional +# dependencies, so the tests that use them are skipped rather than failed when +# they are missing, with a warning so that the gap is not silent. The rest of +# the tests either stub out OpenMM or do not get as far as needing it, and so +# run either way. +OPENMM_SKIP_MSG = 'PDBFixer and OpenMM are needed to fix a structure' +try: + import openmm + import pdbfixer + HAVE_OPENMM = True +except ImportError: + HAVE_OPENMM = False + warnings.warn(OPENMM_SKIP_MSG + ', so the ClustENM tests that need them ' + 'are being skipped') + + +def mockOpenMMApp(): + """Returns a ``patch.dict`` of :mod:`sys.modules` and a mock standing in + for :class:`openmm.app.PDBFile`, so that tests of code which only writes + files can run without OpenMM installed.""" + + pdbfile = MagicMock() + app = types.ModuleType('openmm.app') + app.PDBFile = pdbfile + openmm = types.ModuleType('openmm') + openmm.app = app + patcher = patch.dict(sys.modules, {'openmm': openmm, 'openmm.app': app}) + return patcher, pdbfile + + +class TestClustENM(unittest.TestCase): + """Tests of the arguments of :meth:`.ClustENM.run`, which are checked + before any structure or force field is touched. + + Every one of these but the first raises ValueError, so each is matched + against its message to make sure it is not passing for another reason.""" + + def testParallelWrongType(self): + """Test response to wrong type *parallel* argument.""" + if prody.PY3K: + with self.assertRaises( + TypeError, + msg='ClustENM run with a wrong type parallel failed to ' + 'raise a TypeError'): + ClustENM().run(parallel='nogood') + + def testNoMaxclustOrThreshold(self): + """Test response to neither *maxclust* nor *threshold* being set. + + Only generations that sample conformers need to be clustered, so this + is not required when there are none, as TestClustENMMinimiseOnly + covers.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'maxclust or threshold', + msg='ClustENM run without maxclust or threshold failed to ' + 'say that one of them is needed'): + ClustENM().run() + + def testMaxclustSizeMismatch(self): + """Test response to *maxclust* not covering every generation.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'maxclusts were given', + msg='ClustENM run with too few maxclusts failed to say so'): + ClustENM().run(n_gens=5, maxclust=(10, 30)) + + def testThresholdSizeMismatch(self): + """Test response to *threshold* not covering every generation.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'thresholds were given', + msg='ClustENM run with too few thresholds failed to say so'): + ClustENM().run(n_gens=5, threshold=(1.5, 2.0)) + + def testAtomsUnset(self): + """Test response to running without having set atoms.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'atoms are not set', + msg='ClustENM run without atoms failed to say so'): + ClustENM().run(maxclust=2) + + def testRunWhenBuilt(self): + """Test response to running an ensemble that is already built.""" + if prody.PY3K: + clustenm = ClustENM() + # standing in for a run, which _isBuilt tests for by its conformers + clustenm._confs = np.zeros((1, 3, 3)) + with self.assertRaisesRegex( + ValueError, 'has been built', + msg='ClustENM run on a built ensemble failed to say so'): + clustenm.run(maxclust=2) + + def testEmptyIndex(self): + """Test response to indexing with an empty tuple.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'cannot be empty', + msg='ClustENM indexed with an empty tuple failed to say so'): + ClustENM()[()] + + +class TestClustENMUnbuilt(unittest.TestCase): + """Tests of a ClustENM that has not been given atoms or run yet.""" + + def setUp(self): + if prody.PY3K: + self.clustenm = ClustENM() + + def testNotBuilt(self): + """Test that a new ClustENM is not built""" + if prody.PY3K: + self.assertFalse(self.clustenm._isBuilt(), + 'a new ClustENM claimed to be built') + + def testNoAtoms(self): + """Test that a new ClustENM has no atoms""" + if prody.PY3K: + self.assertIsNone(self.clustenm.getAtoms(), + 'a new ClustENM had atoms') + + def testNoConfs(self): + """Test that a new ClustENM has no conformers""" + if prody.PY3K: + assert_equal(self.clustenm.numConfs(), 0, + 'a new ClustENM had conformers') + + def testNoKeys(self): + """Test that a new ClustENM has no generation keys""" + if prody.PY3K: + self.assertIsNone(self.clustenm.getKeys(), + 'a new ClustENM had generation keys') + + def testDefaultGenerations(self): + """Test the default number of generations""" + if prody.PY3K: + assert_equal(self.clustenm.numGenerations(), 5, + 'ClustENM failed to default to 5 generations') + + +class TestClustENMTitle(unittest.TestCase): + + def testTitleWrongType(self): + """Test response to wrong type *title* argument.""" + if prody.PY3K: + with self.assertRaises( + TypeError, + msg='ClustENM given a wrong type title failed to raise a ' + 'TypeError'): + ClustENM().setTitle(1) + + def testTitleUnset(self): + """Test title of a ClustENM without atoms or an explicit title""" + if prody.PY3K: + assert_equal(ClustENM().getTitle(), 'Unknown', + 'ClustENM without atoms failed to give a dummy title') + + def testTitleRoundTrip(self): + """Test that a title survives being set""" + if prody.PY3K: + clustenm = ClustENM() + clustenm.setTitle(CLUSTENM_TITLE) + assert_equal(clustenm.getTitle(), CLUSTENM_TITLE, + 'ClustENM failed to give back the title it was set') + + +class TestClustENMBlocks(unittest.TestCase): + """The block based subclasses cannot build their ANM without blocks, and + check for them before anything else. That check must not stand in for the + ones their parent makes, so each of those is tested again here with blocks + set, where it would otherwise be masked.""" + + BLOCK_CLASSES = (ClustRTB, ClustImANM) + + def withBlocks(self, cls): + """Returns an instance of *cls* whose blocks are set, so that its own + check passes and those of :meth:`.ClustENM.run` are reached. The blocks + themselves are never used, because no run here gets as far as building + an ANM.""" + + clustenm = cls() + clustenm.setBlocks(np.zeros(10, dtype=int)) + return clustenm + + def testBlocksUnset(self): + """Test response to running without blocks.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'blocks are not set', + msg='%s run without blocks failed to say so' + % cls.__name__): + cls().run(maxclust=2) + + def testBlocksSetReachesAtomsCheck(self): + """Test that blocks being set does not mask unset atoms.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'atoms are not set', + msg='%s run with blocks but no atoms failed to say ' + 'that the atoms are missing' % cls.__name__): + self.withBlocks(cls).run(maxclust=2) + + def testBlocksSetReachesParallelCheck(self): + """Test that blocks being set does not mask a wrong type *parallel*.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaises( + TypeError, + msg='%s run with blocks and a wrong type parallel ' + 'failed to raise a TypeError' % cls.__name__): + self.withBlocks(cls).run(maxclust=2, parallel='nogood') + + def testBlocksSetReachesClusteringCheck(self): + """Test that blocks being set does not mask unset clustering.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'maxclust or threshold', + msg='%s run with blocks but no maxclust or threshold ' + 'failed to say that one of them is needed' + % cls.__name__): + self.withBlocks(cls).run() + + +class TestClustENMBlocksMatchAtoms(unittest.TestCase): + """The blocks of the block based subclasses have to describe the coarse + grained nodes of the structure, whose number is only known once atoms are + set, so these need a really fixed structure to check against. What counts + as valid blocks is checked without one in TestCheckBlocks, and against RTB + and imANM directly in TestRTBBlocks.""" + + BLOCK_CLASSES = (ClustRTB, ClustImANM) + + def setUp(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + self.cwd = os.getcwd() + self.tmpdir = tempfile.mkdtemp() + os.chdir(self.tmpdir) + + self.ATOMS = parseDatafile('1ubi') + + def tearDown(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(self.cwd) + shutil.rmtree(self.tmpdir, ignore_errors=True) + + def runWithBlocks(self, cls, blocks): + """Runs *cls* over the fixed structure with *blocks*, which is given the + number of coarse grained nodes to size itself against.""" + + clustenm = cls() + clustenm.setAtoms(self.ATOMS) + clustenm.setBlocks(blocks(clustenm._n_cg)) + clustenm.run(n_confs=2, n_gens=1, n_modes=2, maxclust=2, + sim=False, outlier=False) + return clustenm + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testBlocksWrongLength(self): + """Test response to blocks that do not cover every node. + + The nodes the blocks are counted against are the coarse grained ones, + not every atom of the fixed structure, so the message says so.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'must match number of coarse grained nodes', + msg='%s run with too few blocks failed to say so' + % cls.__name__): + self.runWithBlocks( + cls, lambda n: np.zeros(n - 5, dtype=int)) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testBlocksAllOneValue(self): + """Test response to blocks that are all the same value. + + One block makes the structure a single rigid body, which has no + internal modes to sample along, so this has to be rejected rather than + left to fail in the eigendecomposition downstream.""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + with self.assertRaisesRegex( + ValueError, 'degrees of freedom of a rigid body', + msg='%s run with a single block failed to say so' + % cls.__name__): + self.runWithBlocks(cls, lambda n: np.zeros(n, dtype=int)) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testBlocksTwoValues(self): + """Test that two blocks are enough to sample along""" + if prody.PY3K: + for cls in self.BLOCK_CLASSES: + clustenm = self.runWithBlocks( + cls, lambda n: np.arange(n) // (n // 2)) + assert_equal(clustenm.numConfs(), 2, + '%s run with two blocks failed to give the ' + 'conformers of both generations' % cls.__name__) + + +class TestClustENMWritePDBFixed(unittest.TestCase): + """Tests for the filename handling of :meth:`.ClustENM.writePDBFixed`, + which :func:`.runANMD` depends on to find the fixed structure it wrote. + A title containing spaces used to make the two disagree, because ANMD + replaced the spaces with underscores and writePDBFixed did not.""" + + def setUp(self): + if not prody.PY3K: + return + + # writePDBFixed writes into the current directory + self.cwd = os.getcwd() + self.tmpdir = tempfile.mkdtemp() + os.chdir(self.tmpdir) + + self.patcher, self.PDBFile = mockOpenMMApp() + self.patcher.start() + + self.clustenm = ClustENM() + self.clustenm.setTitle(CLUSTENM_TITLE) + + def tearDown(self): + if not prody.PY3K: + return + + self.patcher.stop() + os.chdir(self.cwd) + shutil.rmtree(self.tmpdir, ignore_errors=True) + + def testDefaultFilename(self): + """Test that the default filename is derived from the title""" + if prody.PY3K: + filename = self.clustenm.writePDBFixed() + assert_equal(filename, SPACED_NAME, + 'writePDBFixed failed to derive the default filename ' + 'from the title') + self.assertTrue(os.path.exists(SPACED_NAME), + 'writePDBFixed failed to write the file it named') + + def testReturnsFilenameWritten(self): + """Test that the returned filename is the one written to""" + if prody.PY3K: + filename = self.clustenm.writePDBFixed(UNDERSCORED_NAME) + self.assertIsNotNone(filename, + 'writePDBFixed failed to return a filename') + self.assertTrue(os.path.exists(filename), + 'writePDBFixed returned a filename it did not write') + + def testExplicitFilename(self): + """Test that an explicit filename is used verbatim""" + if prody.PY3K: + filename = self.clustenm.writePDBFixed(UNDERSCORED_NAME) + assert_equal(filename, UNDERSCORED_NAME, + 'writePDBFixed failed to use the filename it was given') + self.assertTrue(os.path.exists(UNDERSCORED_NAME), + 'writePDBFixed failed to write to the filename it ' + 'was given') + + def testExplicitFilenameOverridesTitle(self): + """Test that an explicit filename suppresses the title-derived one. + + This is the case that made runANMD raise FileNotFoundError: it asked + for the underscored name and the spaced one was written instead.""" + if prody.PY3K: + self.clustenm.writePDBFixed(UNDERSCORED_NAME) + self.assertFalse(os.path.exists(SPACED_NAME), + 'writePDBFixed wrote the title-derived filename ' + 'instead of the one it was given') + + def testReplaceSpaces(self): + """Test that *replace_spaces* underscores the default filename""" + if prody.PY3K: + filename = self.clustenm.writePDBFixed(replace_spaces=True) + assert_equal(filename, UNDERSCORED_NAME, + 'writePDBFixed with replace_spaces failed to replace ' + 'the spaces of the default filename') + self.assertTrue(os.path.exists(UNDERSCORED_NAME), + 'writePDBFixed with replace_spaces failed to write ' + 'the file it named') + self.assertFalse(os.path.exists(SPACED_NAME), + 'writePDBFixed with replace_spaces wrote a ' + 'filename containing spaces') + + def testReplaceSpacesWithFilename(self): + """Test that *replace_spaces* also applies to an explicit filename""" + if prody.PY3K: + filename = self.clustenm.writePDBFixed(SPACED_NAME, + replace_spaces=True) + assert_equal(filename, UNDERSCORED_NAME, + 'writePDBFixed with replace_spaces failed to replace ' + 'the spaces of the given filename') + + def testWriteFileArguments(self): + """Test that the fixed topology and positions are the ones written""" + if prody.PY3K: + self.clustenm.writePDBFixed(UNDERSCORED_NAME) + assert_equal(self.PDBFile.writeFile.call_count, 1, + 'writePDBFixed failed to write the file once') + + args, kwargs = self.PDBFile.writeFile.call_args + assert_equal(args[0], self.clustenm._topology, + 'writePDBFixed failed to write the fixed topology') + assert_equal(args[1], self.clustenm._positions, + 'writePDBFixed failed to write the fixed positions') + assert_equal(kwargs.get('keepIds'), True, + 'writePDBFixed failed to keep the original ids') + + def testStreamClosed(self): + """Test that the file written to is closed afterwards""" + if prody.PY3K: + self.clustenm.writePDBFixed(UNDERSCORED_NAME) + stream = self.PDBFile.writeFile.call_args[0][2] + self.assertTrue(stream.closed, + 'writePDBFixed failed to close the file it wrote') + + +class TestClustENMFixed(unittest.TestCase): + """Tests against a really fixed structure, rather than the stubbed OpenMM + the tests above use, so that fixing itself and the file it is written to + are both checked.""" + + def setUp(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + # writePDBFixed writes into the current directory + self.cwd = os.getcwd() + self.tmpdir = tempfile.mkdtemp() + os.chdir(self.tmpdir) + + self.ATOMS = parseDatafile('1ubi') + self.ATOMS.setTitle(TITLE) + + def tearDown(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(self.cwd) + shutil.rmtree(self.tmpdir, ignore_errors=True) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testSetAtomsFixes(self): + """Test that setAtoms fixes the structure it is given""" + if prody.PY3K: + clustenm = ClustENM() + clustenm.setAtoms(self.ATOMS) + + assert_equal(clustenm.getTitle(), CLUSTENM_TITLE, + 'setAtoms failed to take the title of the atoms') + # fixing adds the hydrogens that 1ubi does not have + self.assertGreater(clustenm.getAtoms().numAtoms(), + self.ATOMS.numAtoms(), + 'setAtoms failed to add the missing atoms') + self.assertIsNotNone(clustenm.getAtoms().hydrogen, + 'setAtoms failed to add the missing hydrogens') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testWritePDBFixedSpacedTitle(self): + """Test writing a really fixed structure whose title has a space""" + if prody.PY3K: + clustenm = ClustENM() + clustenm.setAtoms(self.ATOMS) + + filename = clustenm.writePDBFixed(UNDERSCORED_NAME) + assert_equal(filename, UNDERSCORED_NAME, + 'writePDBFixed failed to use the filename it was given') + self.assertTrue(os.path.exists(filename), + 'writePDBFixed returned a filename it did not write') + + fixed = parsePDB(filename, compressed=False) + self.assertIsNotNone(fixed, + 'writePDBFixed failed to write a PDB file') + assert_equal(fixed.numAtoms(), clustenm.getAtoms().numAtoms(), + 'writePDBFixed failed to write every fixed atom') + + +class TestClustENMRun(unittest.TestCase): + """Tests of a small ClustENM run. Sampling and clustering are cheap, so + *sim* is turned off to skip the molecular dynamics, leaving the + minimisation of each conformer as the only expensive part.""" + + N_CONFS = 2 + N_GENS = 1 + + @classmethod + def setUpClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + cls.cwd = os.getcwd() + cls.tmpdir = tempfile.mkdtemp() + os.chdir(cls.tmpdir) + + cls.ATOMS = parseDatafile('1ubi') + cls.CLUSTENM = ClustENM() + cls.CLUSTENM.setAtoms(cls.ATOMS) + cls.CLUSTENM.run(n_confs=cls.N_CONFS, n_gens=cls.N_GENS, n_modes=2, + maxclust=2, sim=False, outlier=False) + + @classmethod + def tearDownClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(cls.cwd) + shutil.rmtree(cls.tmpdir, ignore_errors=True) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testBuilt(self): + """Test that a run builds the ensemble""" + if prody.PY3K: + self.assertTrue(self.CLUSTENM._isBuilt(), + 'ClustENM run failed to build the ensemble') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testGenerations(self): + """Test that a run gives the number of generations asked for""" + if prody.PY3K: + assert_equal(self.CLUSTENM.numGenerations(), self.N_GENS, + 'ClustENM run failed to give the number of ' + 'generations it was asked for') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testNumConfs(self): + """Test the conformers of each generation of a run""" + if prody.PY3K: + # generation 0 is the minimised starting structure + assert_equal(self.CLUSTENM.numConfs(0), 1, + 'ClustENM run failed to give 1 conformer for ' + 'generation 0') + assert_equal(self.CLUSTENM.numConfs(1), self.N_CONFS, + 'ClustENM run failed to give %d conformers for ' + 'generation 1' % self.N_CONFS) + assert_equal(self.CLUSTENM.numConfs(), 1 + self.N_CONFS, + 'ClustENM run failed to give the conformers of every ' + 'generation') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testKeysAndLabels(self): + """Test that every conformer is keyed and labelled by generation""" + if prody.PY3K: + keys = [list(key) for key in self.CLUSTENM.getKeys()] + assert_equal(keys, [[0, 0], [1, 0], [1, 1]], + 'ClustENM run failed to key the conformers by ' + 'generation') + assert_equal(self.CLUSTENM.getLabels(), ['0_0', '1_0', '1_1'], + 'ClustENM run failed to label the conformers by ' + 'generation') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testCoordsets(self): + """Test that a run gives coordinates for every fixed atom""" + if prody.PY3K: + coordsets = self.CLUSTENM.getCoordsets() + assert_equal(coordsets.shape, + (1 + self.N_CONFS, + self.CLUSTENM.getAtoms().numAtoms(), 3), + 'ClustENM run failed to give coordinates for every ' + 'conformer and fixed atom') + self.assertTrue(np.isfinite(coordsets).all(), + 'ClustENM run gave coordinates that are not finite') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testPotentials(self): + """Test that a run gives a potential energy for every conformer""" + if prody.PY3K: + potentials = self.CLUSTENM.getPotentials() + assert_equal(len(potentials), 1 + self.N_CONFS, + 'ClustENM run failed to give a potential energy for ' + 'every conformer') + self.assertTrue(np.isfinite(potentials).all(), + 'ClustENM run gave potential energies that are ' + 'not finite') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testRerun(self): + """Test response to running an ensemble that has been run already""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'has been built', + msg='ClustENM run again after a run failed to say that it ' + 'has been built already'): + self.CLUSTENM.run(maxclust=2) + + +class TestClustENMMinimiseOnly(unittest.TestCase): + """Zero generations is a supported way to minimise the starting structure + without generating any conformers along combinations of normal modes. No + generation is clustered, so *maxclust* and *threshold* are not needed + either, which is why this run gives neither.""" + + @classmethod + def setUpClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + cls.cwd = os.getcwd() + cls.tmpdir = tempfile.mkdtemp() + os.chdir(cls.tmpdir) + + cls.ATOMS = parseDatafile('1ubi') + cls.CLUSTENM = ClustENM() + cls.CLUSTENM.setAtoms(cls.ATOMS) + cls.CLUSTENM.run(n_gens=0, n_modes=2, sim=False, outlier=False) + + @classmethod + def tearDownClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(cls.cwd) + shutil.rmtree(cls.tmpdir, ignore_errors=True) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testBuilt(self): + """Test that a run without generations still builds the ensemble""" + if prody.PY3K: + self.assertTrue(self.CLUSTENM._isBuilt(), + 'ClustENM run with 0 generations failed to build ' + 'the ensemble') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testNoGenerations(self): + """Test that no generations are sampled""" + if prody.PY3K: + assert_equal(self.CLUSTENM.numGenerations(), 0, + 'ClustENM run with 0 generations failed to give 0 ' + 'generations') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testOnlyStartingStructure(self): + """Test that only the minimised starting structure is given""" + if prody.PY3K: + assert_equal(self.CLUSTENM.numConfs(), 1, + 'ClustENM run with 0 generations failed to give only ' + 'the minimised starting structure') + + keys = [list(key) for key in self.CLUSTENM.getKeys()] + assert_equal(keys, [[0, 0]], + 'ClustENM run with 0 generations failed to key the ' + 'starting structure as generation 0') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testMinimised(self): + """Test that the starting structure was minimised""" + if prody.PY3K: + coordsets = self.CLUSTENM.getCoordsets() + assert_equal(coordsets.shape, + (1, self.CLUSTENM.getAtoms().numAtoms(), 3), + 'ClustENM run with 0 generations failed to give ' + 'coordinates for every fixed atom') + self.assertTrue(np.isfinite(coordsets).all(), + 'ClustENM run with 0 generations gave coordinates ' + 'that are not finite') + + potentials = self.CLUSTENM.getPotentials() + assert_equal(len(potentials), 1, + 'ClustENM run with 0 generations failed to give a ' + 'potential energy for the starting structure') + self.assertTrue(np.isfinite(potentials).all(), + 'ClustENM run with 0 generations gave a potential ' + 'energy that is not finite') + + +if __name__ == '__main__': + unittest.main() From 500c968db1fe2fe9b383c31d77905b2d1e07b7fd Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 16:46:03 +0100 Subject: [PATCH 04/11] fix: import OpenMM from its current namespace in runANMD runANMD imported OpenMM only through the legacy simtk namespace, which is a deprecated shim, and did so before validating its arguments. So without OpenMM installed the five argument type tests failed with ImportError: Please install PDBFixer and OpenMM to use ANMD instead of the TypeError they assert, because rejecting a bad argument had been made to depend on a simulation package it does not need. The arguments are now validated first, and the import tries the modern openmm namespace before falling back to simtk, as addMissingAtoms already did. ANMD was the last place importing only the legacy one. This also drops the deprecation warning the shim emits. OpenMM and PDBFixer are optional, so the tests that run ANMD are now skipped rather than failed when they are missing, with a warning so the gap is not silent. The argument tests need neither and so always run. Co-Authored-By: Claude Opus 5 --- prody/dynamics/anmd.py | 61 ++++++++++++++++++------- prody/tests/dynamics/test_anmd.py | 75 ++++++++++++++++++++++++++++++- 2 files changed, 118 insertions(+), 18 deletions(-) diff --git a/prody/dynamics/anmd.py b/prody/dynamics/anmd.py index 37c12640f..03ac6b547 100755 --- a/prody/dynamics/anmd.py +++ b/prody/dynamics/anmd.py @@ -25,14 +25,11 @@ __credits__ = ['James Krieger'] __email__ = ['anupam.banerjee@stonybrook.edu', 'jamesmkrieger@gmail.com'] -from numbers import Number import os +from numbers import Number from prody import LOGGER from prody.atomic.atomic import Atomic -from prody.ensemble.ensemble import Ensemble -from prody.proteins.pdbfile import parsePDB, writePDB - from prody.dynamics.anm import ANM from prody.dynamics.clustenm import ClustENM from prody.dynamics.editing import extendModel @@ -40,7 +37,8 @@ from prody.dynamics.nma import NMA from prody.dynamics.pca import PCA from prody.dynamics.sampling import traverseMode - +from prody.ensemble.ensemble import Ensemble +from prody.proteins.pdbfile import parsePDB, writePDB __all__ = ['runANMD'] @@ -89,15 +87,6 @@ def runANMD(atoms, num_modes=2, max_rmsd=2., num_steps=5, tolerance=10.0, A molecular assessment of the alterations in the spike-host protein interactions. *iScience* **2022** 25(3):103939. """ - try: - from simtk.openmm.app import PDBFile, ForceField, \ - Simulation, HBonds, NoCutoff - from simtk.openmm import LangevinIntegrator - from simtk.unit import nanometer, kelvin, picosecond, picoseconds, \ - angstrom, kilojoule, mole - except ImportError: - raise ImportError('Please install PDBFixer and OpenMM to use ANMD') - if not isinstance(atoms, Atomic): raise TypeError('atoms should be an Atomic object') @@ -112,7 +101,6 @@ def runANMD(atoms, num_modes=2, max_rmsd=2., num_steps=5, tolerance=10.0, if not isinstance(tolerance, Number): raise TypeError('tolerance should be a float') - tolerance = tolerance * kilojoule/mole/nanometer pos = kwargs.get('pos', True) if not isinstance(pos, bool): @@ -134,6 +122,44 @@ def runANMD(atoms, num_modes=2, max_rmsd=2., num_steps=5, tolerance=10.0, if not isinstance(anm, (type(None), NMA, ModeSet)): raise TypeError('anm should be an NMA or ModeSet object') + # OpenMM 7.6 moved these out of the legacy simtk namespace, which is + # still shipped as a deprecated shim, so try the modern names first. + try: + try: + from openmm import LangevinIntegrator + from openmm.app import ForceField, HBonds, NoCutoff, PDBFile, Simulation + from openmm.unit import ( + angstrom, + kelvin, + kilojoule, + mole, + nanometer, + picosecond, + picoseconds, + ) + except ImportError: + from simtk.openmm import LangevinIntegrator + from simtk.openmm.app import ( + ForceField, + HBonds, + NoCutoff, + PDBFile, + Simulation, + ) + from simtk.unit import ( + angstrom, + kelvin, + kilojoule, + mole, + nanometer, + picosecond, + picoseconds, + ) + except ImportError: + raise ImportError('Please install PDBFixer and OpenMM to use ANMD') + + tolerance = tolerance * kilojoule/mole/nanometer + pdb_name=atoms.getTitle().replace(' ', '_') fix_name = pdb_name + '_fixed.pdb' @@ -142,7 +168,7 @@ def runANMD(atoms, num_modes=2, max_rmsd=2., num_steps=5, tolerance=10.0, else: clustenm=ClustENM() clustenm.setAtoms(atoms) - clustenm.writePDBFixed() + fix_name = clustenm.writePDBFixed(fix_name) pdb_fix = PDBFile(fix_name) fixmin_name=pdb_name + '_fixedmin.pdb' @@ -163,7 +189,8 @@ def runANMD(atoms, num_modes=2, max_rmsd=2., num_steps=5, tolerance=10.0, simulation.context.setPositions(pdb_fix.positions) simulation.minimizeEnergy(tolerance=tolerance) positions = simulation.context.getState(getPositions=True).getPositions() - PDBFile.writeFile(simulation.topology, positions, open(fixmin_name, 'w')) + with open(fixmin_name, 'w') as stream: + PDBFile.writeFile(simulation.topology, positions, stream) LOGGER.report('The fixed structure was minimised in %.2fs.\n', label='_anmd_min') diff --git a/prody/tests/dynamics/test_anmd.py b/prody/tests/dynamics/test_anmd.py index 9f21b4a77..7f76e7a63 100644 --- a/prody/tests/dynamics/test_anmd.py +++ b/prody/tests/dynamics/test_anmd.py @@ -1,10 +1,16 @@ """This module contains unit tests for :mod:`~prody.dynamics`.""" import os +import shutil import sys +import tempfile +import warnings + import numpy as np from numpy.testing import * + from prody.utilities import importDec + dec = importDec() from prody import * @@ -16,7 +22,7 @@ try: from unittest.mock import MagicMock, patch except ImportError: - from mock import MagicMock, patch + from mock import MagicMock, patch # noqa: UP026 # Prevent threading hangs on remote servers os.environ['OMP_NUM_THREADS'] = '1' @@ -29,6 +35,19 @@ ENSEMBLE = PDBEnsemble(parseDatafile('anmd')) ENSEMBLE.setCoords(ENSEMBLE.getCoordsets()[2]) +# runANMD needs both of these to get as far as minimising anything. They are +# optional dependencies, so the tests that need them are skipped rather than +# failed when they are missing, with a warning so that the gap is not silent. +OPENMM_SKIP_MSG = 'PDBFixer and OpenMM are needed to run ANMD' +try: + import openmm + import pdbfixer + HAVE_OPENMM = True +except ImportError: + HAVE_OPENMM = False + warnings.warn(OPENMM_SKIP_MSG + ', so the ANMD tests that run it are ' + 'being skipped') + class TestANMD(unittest.TestCase): def setUp(self): @@ -141,5 +160,59 @@ def testResultsNumModes1(self): rtol=1e-10, atol=0.25, # may not be so close err_msg='runANMD with num_modes=1 failed to give expected RMSDs') +class TestAnmdSpacedTitle(unittest.TestCase): + """runANMD hands the fixed structure to OpenMM through a PDB file whose + name it derives from the title of *atoms*, replacing spaces with + underscores. It used to look for that name while ClustENM wrote the + unreplaced one, so a title containing a space raised FileNotFoundError.""" + + def setUp(self): + if not prody.PY3K: + return + + # runANMD writes its intermediate PDB files into the current directory + self.cwd = os.getcwd() + self.tmpdir = tempfile.mkdtemp() + os.chdir(self.tmpdir) + + self.ATOMS = parseDatafile('1ubi') + self.ATOMS.setTitle('my prot') + + def tearDown(self): + if not prody.PY3K: + return + + os.chdir(self.cwd) + shutil.rmtree(self.tmpdir, ignore_errors=True) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testAnmdSpacedTitle(self): + """Test that a title containing a space is handled""" + if prody.PY3K: + # a loose tolerance keeps the minimisation short + ensembles = runANMD(self.ATOMS, num_modes=1, num_steps=1, + tolerance=100.) + + assert_equal(len(ensembles), 1, + 'runANMD with a spaced title failed to give 1 ensemble') + assert_equal(len(ensembles[0]), 3, + 'runANMD with a spaced title failed to give an ' + 'ensemble with 3 conformers') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testAnmdSpacedTitleCleansUp(self): + """Test that the intermediate files of a spaced title are removed""" + if prody.PY3K: + runANMD(self.ATOMS, num_modes=1, num_steps=1, tolerance=100.) + + leftover = [name for name in os.listdir(self.tmpdir) + if name.endswith('.pdb')] + assert_equal(leftover, [], + 'runANMD with a spaced title left intermediate files ' # noqa: UP031 + 'behind: %s' % leftover) + + if __name__ == '__main__': unittest.main() From 556e175ce2936cf841e9c25a1034f24396acab6e Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 17:23:27 +0100 Subject: [PATCH 05/11] test: clean up the files the tests leave behind The interaction tests removed three of the four .npy files they save, leaving test_3o21_disu.npy in the working directory after every run, and removed them unconditionally, so a test that failed before saving one raised FileNotFoundError from tearDownClass and reported that in place of the failure which caused it. The buildMSA tests removed nothing at all. buildMSA writes the sequences it aligns beside its output, both named after the title, and clustalw adds an .aln alignment and a .dnd guide tree of its own, so a run left Unknown.fasta, Unknown.aln and Unknown.dnd behind. Co-Authored-By: Claude Opus 5 --- prody/tests/proteins/test_insty.py | 8 ++++++-- prody/tests/sequence/test_analysis.py | 13 +++++++++++++ 2 files changed, 19 insertions(+), 2 deletions(-) diff --git a/prody/tests/proteins/test_insty.py b/prody/tests/proteins/test_insty.py index edfa080a4..ece9b95cf 100644 --- a/prody/tests/proteins/test_insty.py +++ b/prody/tests/proteins/test_insty.py @@ -260,5 +260,9 @@ def testImportHpb(self): def tearDownClass(cls): if prody.PY3K: import os - for filename in ['test_2k39_all.npy', 'test_2k39_sbs.npy', 'test_2k39_disu.npy']: - os.remove(filename) + for filename in ['test_2k39_all.npy', 'test_2k39_sbs.npy', + 'test_2k39_disu.npy', 'test_3o21_disu.npy']: + # a test that failed before saving leaves nothing to remove, + # and raising here would report this instead of that failure + if os.path.isfile(filename): + os.remove(filename) diff --git a/prody/tests/sequence/test_analysis.py b/prody/tests/sequence/test_analysis.py index 73ec68222..f7d74f47b 100644 --- a/prody/tests/sequence/test_analysis.py +++ b/prody/tests/sequence/test_analysis.py @@ -1,5 +1,7 @@ __author__ = 'Ahmet Bakan, Anindita Dutta, Wenzhi Mao, James Krieger' +import os + from prody.tests import TestCase from numpy import array, log, zeros, char, ones, fromfile @@ -1201,6 +1203,17 @@ def testMATLAB10(self): class TestBuildMSA(TestCase): + def tearDown(self): + # buildMSA writes the sequences it aligns beside its output, both named + # after the title, and clustalw adds an .aln alignment and a .dnd guide + # tree of its own. A test that failed before writing one leaves nothing + # to remove, and raising here would report that instead of the failure + # which caused it. + for extension in ('.fasta', '.aln', '.dnd'): + filename = 'Unknown' + extension + if os.path.isfile(filename): + os.remove(filename) + def testBuildMSAlocal(self): sequences = [ags[0].protein["A"].getSequence(), ags[1].protein["A"].getSequence()] From 1cc45d4745d785cd14c45e51b1ecf7be4c26b74e Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 17:23:52 +0100 Subject: [PATCH 06/11] test: cover the clustalo, title and outfilename paths of buildMSA buildMSA has aligned with clustalo as well as clustalw for some time, but neither its docstring nor its tests said so, and the same went for the outfilename argument. Documents both and adds clustalo to environment.yml beside clustalw, which bioconda serves for every platform the CI runs on. The new tests cover aligning with clustalo, naming the output after a given title, and writing it to a given outfilename. clustalo gives an alignment of its own rather than the one clustalw gives, 374 columns against 399 for these sequences, so it is checked for being a valid alignment of the sequences given, by its labels, its column count and the sequences it holds once the gaps are taken out, rather than against a stored alignment that would also pin the version. Neither program is a required dependency, so the tests that need them are now skipped rather than failed when they are missing, as the dssp tests already do, with a warning so that the gap is not silent. The Biopython test needs neither and so still always runs. Co-Authored-By: Claude Opus 5 --- environment.yml | 3 +- prody/sequence/analysis.py | 7 +- prody/tests/sequence/test_analysis.py | 112 ++++++++++++++++++++++---- 3 files changed, 103 insertions(+), 19 deletions(-) diff --git a/environment.yml b/environment.yml index 50dd9ab2a..a885246cd 100644 --- a/environment.yml +++ b/environment.yml @@ -20,7 +20,8 @@ dependencies: - mdtraj - openmm - clustalw - + - clustalo + # Pip packages - pip - pip: diff --git a/prody/sequence/analysis.py b/prody/sequence/analysis.py index 19a3c72b6..a4fec99b5 100644 --- a/prody/sequence/analysis.py +++ b/prody/sequence/analysis.py @@ -742,10 +742,15 @@ def buildMSA(sequences, title='Unknown', labels=None, **kwargs): :type align: bool :arg method: alignment method, one of either Biopython 'global', - Biopython 'local', 'clustalw', 'clustalw2', 'clustal' + Biopython 'local', 'clustalw', 'clustalw2', 'clustal', 'clustalo' or another software in your path. Default is 'local' :type align: str + + :arg outfilename: name of the file to write the alignment to, which + clustalw also writes beside as an ``.aln`` file. + Default is *title* with a ``.fasta`` extension + :type outfilename: str """ align = kwargs.get('align', True) diff --git a/prody/tests/sequence/test_analysis.py b/prody/tests/sequence/test_analysis.py index f7d74f47b..59baebd94 100644 --- a/prody/tests/sequence/test_analysis.py +++ b/prody/tests/sequence/test_analysis.py @@ -1,8 +1,9 @@ __author__ = 'Ahmet Bakan, Anindita Dutta, Wenzhi Mao, James Krieger' import os +import warnings -from prody.tests import TestCase +from prody.tests import TestCase, skipUnless from numpy import array, log, zeros, char, ones, fromfile from numpy.testing import assert_array_equal, assert_array_almost_equal @@ -13,9 +14,23 @@ from prody import calcMSAOccupancy, buildSeqidMatrix, uniqueSequences from prody import buildOMESMatrix, buildSCAMatrix, calcMeff, buildMSA from prody import buildDirectInfoMatrix +from prody.utilities import which LOGGER.verbosity = None +# buildMSA can align with clustalw or clustalo, neither of which is a required +# dependency, so the tests that align with them are skipped rather than failed +# when they are missing, with a warning so that the gap is not silent. +CLUSTALW = which('clustalw') or which('clustalw2') +CLUSTALO = which('clustalo') +CLUSTALW_SKIP_MSG = 'clustalw is not found' +CLUSTALO_SKIP_MSG = 'clustalo is not found' +for _program, _found, _message in (('clustalw', CLUSTALW, CLUSTALW_SKIP_MSG), + ('clustalo', CLUSTALO, CLUSTALO_SKIP_MSG)): + if _found is None: + warnings.warn(_message + ', so the buildMSA tests that align with it ' + 'are being skipped') + FASTA = parseMSA(pathDatafile('msa_Cys_knot.fasta')) FASTA_ALPHA = char.isalpha(FASTA._msa) FASTA_UPPER = char.upper(FASTA._msa) @@ -1203,29 +1218,92 @@ def testMATLAB10(self): class TestBuildMSA(TestCase): + LABELS = ["A2", "A3"] + TITLE = 'test_buildmsa_title' + OUTFILENAME = 'test_buildmsa_outfilename.fasta' + + def setUp(self): + self.sequences = [ags[0].protein["A"].getSequence(), + ags[1].protein["A"].getSequence()] + def tearDown(self): # buildMSA writes the sequences it aligns beside its output, both named - # after the title, and clustalw adds an .aln alignment and a .dnd guide - # tree of its own. A test that failed before writing one leaves nothing - # to remove, and raising here would report that instead of the failure - # which caused it. - for extension in ('.fasta', '.aln', '.dnd'): - filename = 'Unknown' + extension - if os.path.isfile(filename): - os.remove(filename) + # after the title unless outfilename says otherwise, and clustalw adds + # an .aln alignment and a .dnd guide tree of its own. A test that failed + # before writing one leaves nothing to remove, and raising here would + # report that instead of the failure which caused it. + for name in ('Unknown', self.TITLE, + os.path.splitext(self.OUTFILENAME)[0]): + for extension in ('.fasta', '.aln', '.dnd'): + filename = name + extension + if os.path.isfile(filename): + os.remove(filename) def testBuildMSAlocal(self): - sequences = [ags[0].protein["A"].getSequence(), - ags[1].protein["A"].getSequence()] - expect1 = parseMSA(pathDatafile('msa_3hsyA_3o21A.fasta')) - result = buildMSA(sequences, method="local", labels=["A2", "A3"]) + result = buildMSA(self.sequences, method="local", labels=self.LABELS) assert result == expect1, "The list of expected buildMSA results did not contain " + result + @skipUnless(CLUSTALW, CLUSTALW_SKIP_MSG) def testBuildMSAclustalw(self): - sequences = [ags[0].protein["A"].getSequence(), - ags[1].protein["A"].getSequence()] - expect1 = parseMSA(pathDatafile('msa_3hsyA_3o21A_clustalw.fasta')) - result = buildMSA(sequences, method="clustalw", labels=["A2", "A3"]) + result = buildMSA(self.sequences, method="clustalw", labels=self.LABELS) assert result == expect1, "The list of expected buildMSA clustalw results did not contain " + result + + @skipUnless(CLUSTALO, CLUSTALO_SKIP_MSG) + def testBuildMSAclustalo(self): + """Test aligning with clustalo, which gives an alignment of its own + rather than the one clustalw gives, so it is checked for being a valid + alignment of the sequences instead of against a stored one.""" + + result = buildMSA(self.sequences, method="clustalo", labels=self.LABELS) + + assert_array_equal(result.getLabels(), self.LABELS, + 'buildMSA with clustalo failed to label the ' + 'sequences it was given') + self.assertEqual(len(result), len(self.sequences), + 'buildMSA with clustalo failed to align every sequence') + for i, sequence in enumerate(self.sequences): + aligned = str(result[i]) + self.assertEqual(len(aligned), result.numResidues(), + 'buildMSA with clustalo gave sequence {0} a ' + 'length other than the alignment'.format(i)) + self.assertEqual(aligned.replace('-', ''), sequence, + 'buildMSA with clustalo changed sequence {0} by ' + 'more than the gaps it inserted'.format(i)) + + @skipUnless(CLUSTALW, CLUSTALW_SKIP_MSG) + def testBuildMSAtitle(self): + """Test that the files written are named after a given title""" + + expect1 = parseMSA(pathDatafile('msa_3hsyA_3o21A_clustalw.fasta')) + result = buildMSA(self.sequences, title=self.TITLE, method="clustalw", + labels=self.LABELS) + + assert result == expect1, "buildMSA with a title did not give the expected alignment" + for extension in ('.fasta', '.aln'): + filename = self.TITLE + extension + self.assertTrue(os.path.isfile(filename), + 'buildMSA failed to name its {0} file after the ' + 'title it was given'.format(extension)) + self.assertFalse(os.path.isfile('Unknown.fasta'), + 'buildMSA with a title still wrote to the name the ' + 'default title gives') + self.assertEqual(parseMSA(self.TITLE + '.fasta'), expect1, + 'buildMSA wrote an alignment to the title named file ' + 'other than the one it returned') + + @skipUnless(CLUSTALW, CLUSTALW_SKIP_MSG) + def testBuildMSAoutfilename(self): + """Test that the alignment is written where outfilename says""" + + expect1 = parseMSA(pathDatafile('msa_3hsyA_3o21A_clustalw.fasta')) + result = buildMSA(self.sequences, method="clustalw", labels=self.LABELS, + outfilename=self.OUTFILENAME) + + assert result == expect1, "buildMSA with outfilename did not give the expected alignment" + self.assertTrue(os.path.isfile(self.OUTFILENAME), + 'buildMSA failed to write the alignment to outfilename') + self.assertEqual(parseMSA(self.OUTFILENAME), expect1, + 'buildMSA wrote an alignment to outfilename other ' + 'than the one it returned') From 3b0d656efe3fe52f5872444d33ba46ba3611bda4 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 17:32:29 +0100 Subject: [PATCH 07/11] fix: only give ClustENM simulations properties of a platform it chose _prep_sim asked for {'Precision': 'single'} when no platform was given, so it passed platform-specific properties without a platform to go with them. OpenMM discarded them silently up to 8.0, but 8.1 onwards raises ValueError: Cannot specify platform-specific properties, because the Platform is not specified which failed every ClustENM run that did not name a platform, and so is seen with an installed OpenMM but not an older one. Nothing to do with the absence of a GPU: the properties asked for single precision, which is what CUDA and OpenCL are given, but they were asked for whatever platform OpenMM would have picked. Properties are now only given for a platform we chose, leaving OpenMM to pick the fastest it has with its own defaults otherwise, which is what it already did in effect. The CUDA, OpenCL and CPU cases are unchanged. Co-Authored-By: Claude Opus 5 --- prody/dynamics/clustenm.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/prody/dynamics/clustenm.py b/prody/dynamics/clustenm.py index 6094fa2b6..cc71d5b8e 100644 --- a/prody/dynamics/clustenm.py +++ b/prody/dynamics/clustenm.py @@ -292,12 +292,14 @@ def _prep_sim(self, coords, external_forces=[]): 0.002*picosecond) # precision could be mixed, but single is okay. - platform = self._platform if self._platform is None else Platform.getPlatformByName(self._platform) + platform = None if self._platform is None else Platform.getPlatformByName(self._platform) + + # properties are platform specific, so they can only be given for a + # platform we chose. Without one OpenMM picks the fastest it has and + # uses its own defaults, and asking for properties as well raises. properties = None - if self._platform is None: - properties = {'Precision': 'single'} - elif self._platform in ['CUDA', 'OpenCL']: + if self._platform in ['CUDA', 'OpenCL']: properties = {'Precision': 'single'} elif self._platform == 'CPU': if self._threads == 0: From cfbdef865dcb34eca0dd6e15efb86ecab6cfa0e2 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 17:39:26 +0100 Subject: [PATCH 08/11] test: run ClustENM on each platform OpenMM offers Covers running without naming a platform, which is the case that regressed, and on the CPU, CUDA and OpenCL platforms, each of which ClustENM gives different properties. Which platforms OpenMM has depends on how it was built and on the hardware it finds, so the tests of the ones it does not have are skipped rather than failed, with a warning for the accelerated ones so that not having tried them is not silent. Reverting the previous commit fails testNoPlatform and leaves the rest passing, so these hold that fix rather than only running beside it. Co-Authored-By: Claude Opus 5 --- prody/tests/dynamics/test_clustenm.py | 100 ++++++++++++++++++++++++++ 1 file changed, 100 insertions(+) diff --git a/prody/tests/dynamics/test_clustenm.py b/prody/tests/dynamics/test_clustenm.py index 237593946..9cfcfafb7 100644 --- a/prody/tests/dynamics/test_clustenm.py +++ b/prody/tests/dynamics/test_clustenm.py @@ -54,6 +54,29 @@ warnings.warn(OPENMM_SKIP_MSG + ', so the ClustENM tests that need them ' 'are being skipped') +# Which platforms OpenMM can run a simulation on depends on how it was built +# and on the hardware it finds, so the tests of each are skipped rather than +# failed when it does not have that one, with a warning for the accelerated +# ones so that not having tried them is not silent. +if HAVE_OPENMM: + from openmm import Platform + PLATFORMS = [Platform.getPlatform(i).getName() + for i in range(Platform.getNumPlatforms())] +else: + PLATFORMS = [] + +GPU_PLATFORMS = ('CUDA', 'OpenCL') +_missing = [name for name in GPU_PLATFORMS if name not in PLATFORMS] +if HAVE_OPENMM and _missing: + warnings.warn('OpenMM has no %s platform, so the ClustENM tests that run ' + 'on it are being skipped' % ' or '.join(_missing)) + + +def platformSkipMsg(platform): + """Returns the reason for skipping a test of *platform*.""" + + return 'OpenMM has no %s platform' % platform + def mockOpenMMApp(): """Returns a ``patch.dict`` of :mod:`sys.modules` and a mock standing in @@ -725,5 +748,82 @@ def testMinimised(self): 'energy that is not finite') +class TestClustENMPlatform(unittest.TestCase): + """Properties can only be given to OpenMM for a platform that was named, + so each platform is run to make sure ClustENM asks for a combination of the + two that OpenMM accepts. Without a platform it picks the fastest it has, + and asking for properties as well used to be silently ignored but now + raises, so the run without one is the case that matters most here. + + These are minimisation only runs, which is enough to build a simulation.""" + + def setUp(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + self.cwd = os.getcwd() + self.tmpdir = tempfile.mkdtemp() + os.chdir(self.tmpdir) + + self.ATOMS = parseDatafile('1ubi') + + def tearDown(self): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(self.cwd) + shutil.rmtree(self.tmpdir, ignore_errors=True) + + def runOnPlatform(self, platform): + """Runs ClustENM over the fixed structure on *platform*, or on whichever + one OpenMM picks if it is **None**, and checks the minimised starting + structure that comes back.""" + + clustenm = ClustENM() + clustenm.setAtoms(self.ATOMS) + + kwargs = {} if platform is None else {'platform': platform} + clustenm.run(n_gens=0, n_modes=2, sim=False, outlier=False, **kwargs) + + named = 'no platform' if platform is None else platform + assert_equal(clustenm.numConfs(), 1, + 'ClustENM run with %s failed to give the minimised ' + 'starting structure' % named) + self.assertTrue(np.isfinite(clustenm.getCoordsets()).all(), + 'ClustENM run with %s gave coordinates that are not ' + 'finite' % named) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testNoPlatform(self): + """Test running without naming a platform""" + if prody.PY3K: + self.runOnPlatform(None) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + @unittest.skipUnless('CPU' in PLATFORMS, platformSkipMsg('CPU')) + def testCPUPlatform(self): + """Test running on the CPU platform, which is given a thread count""" + if prody.PY3K: + self.runOnPlatform('CPU') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + @unittest.skipUnless('CUDA' in PLATFORMS, platformSkipMsg('CUDA')) + def testCUDAPlatform(self): + """Test running on the CUDA platform, which is given a precision""" + if prody.PY3K: + self.runOnPlatform('CUDA') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + @unittest.skipUnless('OpenCL' in PLATFORMS, platformSkipMsg('OpenCL')) + def testOpenCLPlatform(self): + """Test running on the OpenCL platform, which is given a precision""" + if prody.PY3K: + self.runOnPlatform('OpenCL') + + if __name__ == '__main__': unittest.main() From 8d5bd5d5a9ac41cde1b5db018c8b70b64877e1bd Mon Sep 17 00:00:00 2001 From: James Krieger Date: Wed, 12 Aug 2026 17:49:24 +0100 Subject: [PATCH 09/11] test: count the conformers a ClustENM generation keeps, not assume them A generation samples n_confs conformers and then filters and clusters them, so how many it keeps depends on the conformers themselves, and so on the modes the eigensolver of the machine gives. The block based runs keep one of the two here and two on the CI, which failed ClustRTB run with two blocks failed to give the conformers of both generations: ACTUAL: 3 DESIRED: 2 so the count is not a property of the run to be expected in advance. The tests now read it and check what does hold: generation 0 is the minimised starting structure alone, generation 1 keeps at least one conformer and no more than it sampled, and the totals, keys, labels, coordinate sets and potential energies all agree with it. Deriving them keeps the relations between them just as strictly, and still fails if a generation samples nothing or keeps more than it sampled, while not pinning a number that belongs to the machine. Also takes n_confs and maxclust for the block runs from the same attribute the bound is read from, so the two cannot drift apart. Co-Authored-By: Claude Opus 5 --- prody/tests/dynamics/test_clustenm.py | 60 +++++++++++++++++++++------ 1 file changed, 47 insertions(+), 13 deletions(-) diff --git a/prody/tests/dynamics/test_clustenm.py b/prody/tests/dynamics/test_clustenm.py index 9cfcfafb7..4410295ad 100644 --- a/prody/tests/dynamics/test_clustenm.py +++ b/prody/tests/dynamics/test_clustenm.py @@ -298,6 +298,9 @@ class TestClustENMBlocksMatchAtoms(unittest.TestCase): BLOCK_CLASSES = (ClustRTB, ClustImANM) + # how many conformers each generation samples, and the most it can keep + N_CONFS = 2 + def setUp(self): if not prody.PY3K or not HAVE_OPENMM: return @@ -322,8 +325,8 @@ def runWithBlocks(self, cls, blocks): clustenm = cls() clustenm.setAtoms(self.ATOMS) clustenm.setBlocks(blocks(clustenm._n_cg)) - clustenm.run(n_confs=2, n_gens=1, n_modes=2, maxclust=2, - sim=False, outlier=False) + clustenm.run(n_confs=self.N_CONFS, n_gens=1, n_modes=2, + maxclust=self.N_CONFS, sim=False, outlier=False) return clustenm @dec.slow @@ -361,12 +364,29 @@ def testBlocksAllOneValue(self): @dec.slow @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) def testBlocksTwoValues(self): - """Test that two blocks are enough to sample along""" + """Test that two blocks are enough to sample along. + + How many of the conformers sampled from them survive being filtered + and clustered depends on the conformers themselves, and so on the + modes the eigensolver of the machine gives, so this counts them rather + than expecting a number of its own.""" if prody.PY3K: for cls in self.BLOCK_CLASSES: clustenm = self.runWithBlocks( cls, lambda n: np.arange(n) // (n // 2)) - assert_equal(clustenm.numConfs(), 2, + + sampled = clustenm.numConfs(1) + assert_equal(clustenm.numConfs(0), 1, + '%s run with two blocks failed to give the ' + 'minimised starting structure' % cls.__name__) + self.assertGreaterEqual(sampled, 1, + '%s run with two blocks failed to ' + 'sample along them' % cls.__name__) + self.assertLessEqual(sampled, self.N_CONFS, + '%s run with two blocks gave more ' + 'conformers than it sampled' + % cls.__name__) + assert_equal(clustenm.numConfs(), 1 + sampled, '%s run with two blocks failed to give the ' 'conformers of both generations' % cls.__name__) @@ -600,16 +620,25 @@ def testGenerations(self): @dec.slow @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) def testNumConfs(self): - """Test the conformers of each generation of a run""" + """Test the conformers of each generation of a run. + + How many of the conformers sampled survive being filtered and + clustered depends on the conformers themselves, and so on the modes + the eigensolver of the machine gives, so generation 1 is counted + rather than expected to keep every one of them.""" if prody.PY3K: # generation 0 is the minimised starting structure assert_equal(self.CLUSTENM.numConfs(0), 1, 'ClustENM run failed to give 1 conformer for ' 'generation 0') - assert_equal(self.CLUSTENM.numConfs(1), self.N_CONFS, - 'ClustENM run failed to give %d conformers for ' - 'generation 1' % self.N_CONFS) - assert_equal(self.CLUSTENM.numConfs(), 1 + self.N_CONFS, + sampled = self.CLUSTENM.numConfs(1) + self.assertGreaterEqual(sampled, 1, + 'ClustENM run failed to give any ' + 'conformers for generation 1') + self.assertLessEqual(sampled, self.N_CONFS, + 'ClustENM run gave more conformers for ' + 'generation 1 than it sampled') + assert_equal(self.CLUSTENM.numConfs(), 1 + sampled, 'ClustENM run failed to give the conformers of every ' 'generation') @@ -618,11 +647,16 @@ def testNumConfs(self): def testKeysAndLabels(self): """Test that every conformer is keyed and labelled by generation""" if prody.PY3K: + # the starting structure, then however many generation 1 kept + expect = [[0, 0]] + [[1, i] + for i in range(self.CLUSTENM.numConfs(1))] + keys = [list(key) for key in self.CLUSTENM.getKeys()] - assert_equal(keys, [[0, 0], [1, 0], [1, 1]], + assert_equal(keys, expect, 'ClustENM run failed to key the conformers by ' 'generation') - assert_equal(self.CLUSTENM.getLabels(), ['0_0', '1_0', '1_1'], + assert_equal(self.CLUSTENM.getLabels(), + ['%d_%d' % tuple(key) for key in expect], 'ClustENM run failed to label the conformers by ' 'generation') @@ -633,7 +667,7 @@ def testCoordsets(self): if prody.PY3K: coordsets = self.CLUSTENM.getCoordsets() assert_equal(coordsets.shape, - (1 + self.N_CONFS, + (self.CLUSTENM.numConfs(), self.CLUSTENM.getAtoms().numAtoms(), 3), 'ClustENM run failed to give coordinates for every ' 'conformer and fixed atom') @@ -646,7 +680,7 @@ def testPotentials(self): """Test that a run gives a potential energy for every conformer""" if prody.PY3K: potentials = self.CLUSTENM.getPotentials() - assert_equal(len(potentials), 1 + self.N_CONFS, + assert_equal(len(potentials), self.CLUSTENM.numConfs(), 'ClustENM run failed to give a potential energy for ' 'every conformer') self.assertTrue(np.isfinite(potentials).all(), From 91c4e51faa35467c2837d9bd25e97ab45421f7a2 Mon Sep 17 00:00:00 2001 From: James Krieger Date: Thu, 13 Aug 2026 10:31:38 +0100 Subject: [PATCH 10/11] test: skip a platform OpenMM lists but cannot make a context on The macOS runners offer the OpenCL platform without a device it can use, so testOpenCLPlatform was not skipped and failed in the run itself: openmm.OpenMMException: No compatible OpenCL platform is available Being one of the platforms OpenMM lists only says the plugin loaded, so each is now tried, by making a context of a single particle on it, rather than looked up in the list. The warning and the skip reasons say that OpenMM cannot run on it here rather than that it does not have it, since the two are different and only the second is about how it was built. Co-Authored-By: Claude Opus 5 --- prody/tests/dynamics/test_clustenm.py | 60 +++++++++++++++++++-------- 1 file changed, 42 insertions(+), 18 deletions(-) diff --git a/prody/tests/dynamics/test_clustenm.py b/prody/tests/dynamics/test_clustenm.py index 4410295ad..194047fa5 100644 --- a/prody/tests/dynamics/test_clustenm.py +++ b/prody/tests/dynamics/test_clustenm.py @@ -56,26 +56,50 @@ # Which platforms OpenMM can run a simulation on depends on how it was built # and on the hardware it finds, so the tests of each are skipped rather than -# failed when it does not have that one, with a warning for the accelerated -# ones so that not having tried them is not silent. -if HAVE_OPENMM: - from openmm import Platform - PLATFORMS = [Platform.getPlatform(i).getName() - for i in range(Platform.getNumPlatforms())] -else: - PLATFORMS = [] - -GPU_PLATFORMS = ('CUDA', 'OpenCL') -_missing = [name for name in GPU_PLATFORMS if name not in PLATFORMS] -if HAVE_OPENMM and _missing: - warnings.warn('OpenMM has no %s platform, so the ClustENM tests that run ' - 'on it are being skipped' % ' or '.join(_missing)) +# failed when it cannot use that one, with a warning so that not having tried +# them is not silent. +def usablePlatform(name): + """Returns whether OpenMM can make a context on the platform called *name*. + + Being one of the platforms it lists is not enough. The OpenCL plugin loads + on machines that have no device it can use, where the platform is offered + but making a context on it raises, so each one is tried rather than looked + up.""" + + if not HAVE_OPENMM: + return False + + from openmm import Context, Platform, System, VerletIntegrator + + listed = [Platform.getPlatform(i).getName() + for i in range(Platform.getNumPlatforms())] + if name not in listed: + return False + + # a single particle is enough to make a context and find out + system = System() + system.addParticle(1.) + try: + Context(system, VerletIntegrator(0.001), + Platform.getPlatformByName(name)) + except Exception: + return False + + return True + + +PLATFORMS = {name: usablePlatform(name) for name in ('CPU', 'CUDA', 'OpenCL')} +_unusable = sorted(name for name, usable in PLATFORMS.items() if not usable) +if HAVE_OPENMM and _unusable: + warnings.warn('OpenMM cannot run on the %s platform here, so the ClustENM ' + 'tests that run on it are being skipped' + % ' or '.join(_unusable)) def platformSkipMsg(platform): """Returns the reason for skipping a test of *platform*.""" - return 'OpenMM has no %s platform' % platform + return 'OpenMM cannot run on the %s platform here' % platform def mockOpenMMApp(): @@ -836,7 +860,7 @@ def testNoPlatform(self): @dec.slow @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) - @unittest.skipUnless('CPU' in PLATFORMS, platformSkipMsg('CPU')) + @unittest.skipUnless(PLATFORMS['CPU'], platformSkipMsg('CPU')) def testCPUPlatform(self): """Test running on the CPU platform, which is given a thread count""" if prody.PY3K: @@ -844,7 +868,7 @@ def testCPUPlatform(self): @dec.slow @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) - @unittest.skipUnless('CUDA' in PLATFORMS, platformSkipMsg('CUDA')) + @unittest.skipUnless(PLATFORMS['CUDA'], platformSkipMsg('CUDA')) def testCUDAPlatform(self): """Test running on the CUDA platform, which is given a precision""" if prody.PY3K: @@ -852,7 +876,7 @@ def testCUDAPlatform(self): @dec.slow @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) - @unittest.skipUnless('OpenCL' in PLATFORMS, platformSkipMsg('OpenCL')) + @unittest.skipUnless(PLATFORMS['OpenCL'], platformSkipMsg('OpenCL')) def testOpenCLPlatform(self): """Test running on the OpenCL platform, which is given a precision""" if prody.PY3K: From 48db0949df2390ccbb582a86947b59a02c64395b Mon Sep 17 00:00:00 2001 From: James Krieger Date: Thu, 13 Aug 2026 12:16:48 +0100 Subject: [PATCH 11/11] test: cover starting ClustENM from more than one conformer setAtoms now takes a list of structures, which become the generation 0 population rather than the single minimised starting structure one structure gives, so generation 0 holds as many conformers as it was given and the keys and labels run 0_0, 0_1, ... before generation 1 begins. Nothing covered that. Covers the initial population, its keys and labels, the generation sampled from it, starting from three structures rather than two, a list of one behaving as a single structure does, an empty list, and structures that are not the same molecule, whose coordinates cannot fit the topology built from the first. How many conformers the sampled generation keeps is counted rather than expected, as the single start tests do. The models of a local NMR structure stand in for the conformers a run would start from, so nothing is fetched. The tests that only read a run from every model share one, as the single start run tests do. Co-Authored-By: Claude Opus 5 --- prody/tests/dynamics/test_clustenm.py | 162 ++++++++++++++++++++++++++ 1 file changed, 162 insertions(+) diff --git a/prody/tests/dynamics/test_clustenm.py b/prody/tests/dynamics/test_clustenm.py index 194047fa5..0f742d7f3 100644 --- a/prody/tests/dynamics/test_clustenm.py +++ b/prody/tests/dynamics/test_clustenm.py @@ -806,6 +806,168 @@ def testMinimised(self): 'energy that is not finite') +class TestClustENMMultiStart(unittest.TestCase): + """setAtoms takes a list of structures to start a run from more than one + conformer. They become the generation 0 population, each of them minimised, + where a single structure gives a generation 0 of just itself. The topology + is built from the first, so the rest have to be the same molecule. + + The models of an NMR structure stand in for the conformers a run would + start from, so that nothing has to be fetched.""" + + N_CONFS = 2 + N_MODELS = 2 + + @classmethod + def setUpClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + cls.cwd = os.getcwd() + cls.tmpdir = tempfile.mkdtemp() + os.chdir(cls.tmpdir) + + cls.path = pathDatafile('multi_model_truncated') + cls.MODELS = [cls.model(i) for i in range(1, cls.N_MODELS + 1)] + + # the tests which only read a run from every model share this one + cls.MULTISTART = cls.runFrom(cls.MODELS) + + @classmethod + def tearDownClass(cls): + if not prody.PY3K or not HAVE_OPENMM: + return + + os.chdir(cls.cwd) + shutil.rmtree(cls.tmpdir, ignore_errors=True) + + @classmethod + def model(cls, number): + """Returns model *number* of the structure the models come from.""" + + return parsePDB(cls.path, model=number).select('protein and chain A') + + @classmethod + def runFrom(cls, models): + """Runs ClustENM started from *models* and returns it.""" + + clustenm = ClustENM('multistart') + clustenm.setAtoms(models) + clustenm.run(n_gens=1, n_modes=2, n_confs=cls.N_CONFS, rmsd=1.0, + maxclust=cls.N_CONFS, sim=False, outlier=False) + return clustenm + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testInitialStructures(self): + """Test that every structure given starts the run""" + if prody.PY3K: + clustenm = self.MULTISTART + + assert_equal(clustenm.numConfs(0), self.N_MODELS, + 'ClustENM run from %d structures failed to start ' + 'generation 0 from all of them' % self.N_MODELS) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testKeysAndLabels(self): + """Test that the initial structures are keyed as generation 0""" + if prody.PY3K: + clustenm = self.MULTISTART + + # every initial structure, then however many generation 1 kept + expect = ([[0, i] for i in range(self.N_MODELS)] + + [[1, j] for j in range(clustenm.numConfs(1))]) + + keys = [list(key) for key in clustenm.getKeys()] + assert_equal(keys, expect, + 'ClustENM run from %d structures failed to key them ' + 'as generation 0' % self.N_MODELS) + assert_equal(clustenm.getLabels(), + ['%d_%d' % tuple(key) for key in expect], + 'ClustENM run from %d structures failed to label them ' + 'as generation 0' % self.N_MODELS) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testSampledGeneration(self): + """Test that a generation is sampled from the structures given. + + How many conformers it keeps depends on the machine, so this counts + them, as the single start tests do.""" + if prody.PY3K: + clustenm = self.MULTISTART + + sampled = clustenm.numConfs(1) + self.assertGreaterEqual(sampled, 1, + 'ClustENM run from several structures ' + 'failed to sample a generation from them') + self.assertLessEqual(sampled, self.N_CONFS, + 'ClustENM run from several structures gave ' + 'more conformers than it sampled') + assert_equal(clustenm.numConfs(), self.N_MODELS + sampled, + 'ClustENM run from several structures failed to give ' + 'the initial structures and the sampled conformers') + + coordsets = clustenm.getCoordsets() + assert_equal(coordsets.shape, + (clustenm.numConfs(), + clustenm.getAtoms().numAtoms(), 3), + 'ClustENM run from several structures failed to give ' + 'coordinates for every conformer and atom') + self.assertTrue(np.isfinite(coordsets).all(), + 'ClustENM run from several structures gave ' + 'coordinates that are not finite') + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testMoreThanTwoStructures(self): + """Test starting from more structures than two""" + if prody.PY3K: + models = self.MODELS + [self.model(self.N_MODELS + 1)] + clustenm = self.runFrom(models) + + assert_equal(clustenm.numConfs(0), len(models), + 'ClustENM run from %d structures failed to start ' + 'generation 0 from all of them' % len(models)) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testOneStructureInAList(self): + """Test that a list of one structure starts a single conformer run""" + if prody.PY3K: + clustenm = self.runFrom(self.MODELS[:1]) + + assert_equal(clustenm.numConfs(0), 1, + 'ClustENM run from a list of one structure failed to ' + 'give the one minimised starting structure') + + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testNoStructures(self): + """Test response to an empty list of structures.""" + if prody.PY3K: + with self.assertRaisesRegex( + ValueError, 'empty list of structures', + msg='ClustENM given no structures to start from failed to ' + 'say so'): + ClustENM().setAtoms([]) + + @dec.slow + @unittest.skipUnless(HAVE_OPENMM, OPENMM_SKIP_MSG) + def testDifferentMolecules(self): + """Test response to structures that are not the same molecule. + + The topology comes from the first, so the coordinates of the rest have + to fit it.""" + if prody.PY3K: + other = parseDatafile('1ubi').select('protein') + with self.assertRaisesRegex( + ValueError, 'must be the same molecule', + msg='ClustENM given structures of different molecules to ' + 'start from failed to say so'): + ClustENM().setAtoms([self.MODELS[0], other]) + + class TestClustENMPlatform(unittest.TestCase): """Properties can only be given to OpenMM for a platform that was named, so each platform is run to make sure ClustENM asks for a combination of the