From d5e83dccb9931d6c71d14ca3612135a048daa39d Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Thu, 6 Aug 2026 23:09:43 +0530 Subject: [PATCH 1/7] Rotate velocities and forces besides positions during trajectory rotations --- package/CHANGELOG | 7 ++++ package/MDAnalysis/core/groups.py | 9 +++- package/MDAnalysis/transformations/rotate.py | 4 ++ .../MDAnalysisTests/analysis/test_align.py | 41 +++++++++++++++++++ .../MDAnalysisTests/core/test_atomgroup.py | 26 ++++++++++++ .../transformations/test_rotate.py | 33 +++++++++++++++ 6 files changed, 119 insertions(+), 1 deletion(-) diff --git a/package/CHANGELOG b/package/CHANGELOG index 71411414874..24161615fea 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -23,6 +23,13 @@ The rules for this file: * 2.11.0 Fixes + * Fixes + * `AtomGroup.rotate()` and the `rotateby` trajectory transformation now + also rotate velocities and forces besides positions. This also affects + `MDAnalysis.analysis.align.alignto()` and `AlignTraj`, since they + apply their fit via `AtomGroup.rotate()`. Note: this changes existing + behavior, as velocities/forces were previously left untouched by + rotation. (Issue #5421, PR #) * Fix FileLock tests for XTC and TRR: lock file is no longer removed (#5382) * InterRDF now correctly returns bins in parallel (PR #5344) * `Merge()` no longer raises a TypeError on Universes that have a `cmaps` diff --git a/package/MDAnalysis/core/groups.py b/package/MDAnalysis/core/groups.py index f6e03ee9798..c5fc2211101 100644 --- a/package/MDAnalysis/core/groups.py +++ b/package/MDAnalysis/core/groups.py @@ -1610,11 +1610,18 @@ def rotate(self, R, point=(0, 0, 0)): require_translation = bool(np.count_nonzero(point)) if require_translation: atomgroup.translate(-point) - x = atomgroup.universe.trajectory.ts.positions + ts = atomgroup.universe.trajectory.ts + x = ts.positions idx = atomgroup.indices x[idx] = np.dot(x[idx], R.T) if require_translation: atomgroup.translate(point) + if ts.has_velocities: + v = ts.velocities + v[idx] = np.dot(v[idx], R.T) + if ts.has_forces: + f = ts.forces + f[idx] = np.dot(f[idx], R.T) return self diff --git a/package/MDAnalysis/transformations/rotate.py b/package/MDAnalysis/transformations/rotate.py index 4d8fa71d0b1..c23152cd790 100644 --- a/package/MDAnalysis/transformations/rotate.py +++ b/package/MDAnalysis/transformations/rotate.py @@ -200,4 +200,8 @@ def _transform(self, ts): translation = matrix[:3, 3] ts.positions = np.dot(ts.positions, rotation) ts.positions += translation + if ts.has_velocities: + ts.velocities = np.dot(ts.velocities, rotation) + if ts.has_forces: + ts.forces = np.dot(ts.forces, rotation) return ts diff --git a/testsuite/MDAnalysisTests/analysis/test_align.py b/testsuite/MDAnalysisTests/analysis/test_align.py index fbda36b1580..a97c6ba241d 100644 --- a/testsuite/MDAnalysisTests/analysis/test_align.py +++ b/testsuite/MDAnalysisTests/analysis/test_align.py @@ -31,6 +31,7 @@ import pytest from MDAnalysis import SelectionError, SelectionWarning from MDAnalysisTests import executable_not_found +from MDAnalysis.lib import transformations from MDAnalysisTests.datafiles import ( PSF, DCD, @@ -821,3 +822,43 @@ def test_alignto_reorder_atomgroups(): ref = u.atoms[[3, 2, 1, 0]] rmsd = align.alignto(mobile, ref, select="bynum 1-4") assert_allclose(rmsd, (0.0, 0.0)) + + +def test_alignto_rotates_velocities_and_forces(): + mobile = mda.Universe.empty( + 4, trajectory=True, velocities=True, forces=True + ) + reference = mda.Universe.empty(4, trajectory=True) + + mobile.add_TopologyAttr("masses", [1.0, 1.0, 1.0, 1.0]) + reference.add_TopologyAttr("masses", [1.0, 1.0, 1.0, 1.0]) + + ref_pos = np.array( + [[1, 0, 0], [0, 1, 0], [-1, 0, 0], [0, -1, 0]], dtype=np.float64 + ) + reference.atoms.positions = ref_pos + + angle = np.pi / 3 + known_R = transformations.rotation_matrix(angle, [0, 0, 1])[:3, :3] + mobile.atoms.positions = np.dot(ref_pos, known_R.T) + + rng = np.random.RandomState(0) + orig_v = rng.random((4, 3)) + orig_f = rng.random((4, 3)) + mobile.atoms.velocities = orig_v.copy() + mobile.atoms.forces = orig_f.copy() + + mobile_centered = ( + mobile.atoms.positions - mobile.atoms.center_of_geometry() + ) + ref_centered = ref_pos - ref_pos.mean(axis=0) + expected_R, _ = align.rotation_matrix(mobile_centered, ref_centered) + + align.alignto(mobile, reference) + + assert_allclose( + mobile.atoms.velocities, np.dot(orig_v, expected_R.T), atol=1e-6 + ) + assert_allclose( + mobile.atoms.forces, np.dot(orig_f, expected_R.T), atol=1e-6 + ) diff --git a/testsuite/MDAnalysisTests/core/test_atomgroup.py b/testsuite/MDAnalysisTests/core/test_atomgroup.py index 384b14a5baf..79e4efcfe45 100644 --- a/testsuite/MDAnalysisTests/core/test_atomgroup.py +++ b/testsuite/MDAnalysisTests/core/test_atomgroup.py @@ -417,6 +417,32 @@ def test_rotateby(self, u, coords): [-2 * np.cos(angle) + 1, -2 * np.sin(angle), 0], ) + def test_rotate_velocities_forces(self): + u = mda.Universe.empty( + 2, trajectory=True, velocities=True, forces=True + ) + u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) + u.atoms.velocities = np.array([[1, 0, 0], [0, 1, 0]]) + u.atoms.forces = np.array([[0, 0, 1], [1, 1, 0]]) + + orig_v = u.atoms.velocities.copy() + orig_f = u.atoms.forces.copy() + + axis = np.array([0, 0, 1]) + for angle in np.linspace(0, np.pi): + R = transformations.rotation_matrix(angle, axis)[:3, :3] + u.atoms.velocities = orig_v.copy() + u.atoms.forces = orig_f.copy() + u.atoms.rotate(R) + assert_almost_equal(u.atoms.velocities, np.dot(orig_v, R.T)) + assert_almost_equal(u.atoms.forces, np.dot(orig_f, R.T)) + + def test_rotate_no_velocities_forces_does_not_raise(self): + u = mda.Universe.empty(2, trajectory=True) + u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) + R = transformations.rotation_matrix(1, [0, 0, 1])[:3, :3] + u.atoms.rotate(R) + def test_transform_rotation_only(self, u, coords): R = np.eye(3) u.atoms.rotate(R) diff --git a/testsuite/MDAnalysisTests/transformations/test_rotate.py b/testsuite/MDAnalysisTests/transformations/test_rotate.py index 4f8fd9867b5..c003d6535b1 100644 --- a/testsuite/MDAnalysisTests/transformations/test_rotate.py +++ b/testsuite/MDAnalysisTests/transformations/test_rotate.py @@ -182,6 +182,39 @@ def test_rotateby_atomgroup_com_pbc(rotate_universes): assert_array_almost_equal(transformed.positions, ref.positions, decimal=6) +def test_rotateby_velocities_forces(): + u = mda.Universe.empty(2, trajectory=True, velocities=True, forces=True) + u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) + u.atoms.velocities = np.array([[1, 0, 0], [0, 1, 0]]) + u.atoms.forces = np.array([[0, 0, 1], [1, 1, 0]]) + ts = u.trajectory.ts + + orig_v = ts.velocities.copy() + orig_f = ts.forces.copy() + + axis = [0, 0, 1] + point = [0, 0, 0] + angle = 90 + matrix = rotation_matrix(np.deg2rad(angle), axis, point) + rotation = matrix[:3, :3].T + + transformed_ts = rotateby(angle, axis, point=point)(ts) + + assert_array_almost_equal( + transformed_ts.velocities, np.dot(orig_v, rotation), decimal=6 + ) + assert_array_almost_equal( + transformed_ts.forces, np.dot(orig_f, rotation), decimal=6 + ) + + +def test_rotateby_no_velocities_forces_does_not_raise(): + u = mda.Universe.empty(2, trajectory=True) + u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) + ts = u.trajectory.ts + rotateby(90, [0, 0, 1], point=[0, 0, 0])(ts) + + @pytest.mark.parametrize( "ag", ( From df93377c4154786f533c8c263b793e9a39181d03 Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Thu, 6 Aug 2026 23:21:15 +0530 Subject: [PATCH 2/7] update CHANGELOG --- package/CHANGELOG | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/package/CHANGELOG b/package/CHANGELOG index 24161615fea..3e7c3d309d3 100644 --- a/package/CHANGELOG +++ b/package/CHANGELOG @@ -23,13 +23,12 @@ The rules for this file: * 2.11.0 Fixes - * Fixes * `AtomGroup.rotate()` and the `rotateby` trajectory transformation now also rotate velocities and forces besides positions. This also affects `MDAnalysis.analysis.align.alignto()` and `AlignTraj`, since they apply their fit via `AtomGroup.rotate()`. Note: this changes existing behavior, as velocities/forces were previously left untouched by - rotation. (Issue #5421, PR #) + rotation. (Issue #5421, PR #5452) * Fix FileLock tests for XTC and TRR: lock file is no longer removed (#5382) * InterRDF now correctly returns bins in parallel (PR #5344) * `Merge()` no longer raises a TypeError on Universes that have a `cmaps` From 00f4db21fae88d0d8249448935f70f189f16e20f Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Sun, 9 Aug 2026 22:36:04 +0530 Subject: [PATCH 3/7] update docs and tests --- package/MDAnalysis/core/groups.py | 9 ++++++--- .../MDAnalysisTests/analysis/test_align.py | 6 ++++-- .../MDAnalysisTests/core/test_atomgroup.py | 6 ++++-- .../transformations/test_rotate.py | 17 +++++++++++++---- 4 files changed, 27 insertions(+), 11 deletions(-) diff --git a/package/MDAnalysis/core/groups.py b/package/MDAnalysis/core/groups.py index c5fc2211101..c30ede06833 100644 --- a/package/MDAnalysis/core/groups.py +++ b/package/MDAnalysis/core/groups.py @@ -1601,6 +1601,8 @@ def rotate(self, R, point=(0, 0, 0)): rotateby : rotate around given axis and angle MDAnalysis.lib.transformations : module of all coordinate transforms + .. versionchanged: 2.11.0 + Also rotate velocities and forces if present in Timestep. """ R = np.asarray(R) point = np.asarray(point) @@ -1611,17 +1613,18 @@ def rotate(self, R, point=(0, 0, 0)): if require_translation: atomgroup.translate(-point) ts = atomgroup.universe.trajectory.ts + R_T = R.T x = ts.positions idx = atomgroup.indices - x[idx] = np.dot(x[idx], R.T) + x[idx] = np.dot(x[idx], R_T) if require_translation: atomgroup.translate(point) if ts.has_velocities: v = ts.velocities - v[idx] = np.dot(v[idx], R.T) + v[idx] = np.dot(v[idx], R_T) if ts.has_forces: f = ts.forces - f[idx] = np.dot(f[idx], R.T) + f[idx] = np.dot(f[idx], R_T) return self diff --git a/testsuite/MDAnalysisTests/analysis/test_align.py b/testsuite/MDAnalysisTests/analysis/test_align.py index a97c6ba241d..5ac553c16c5 100644 --- a/testsuite/MDAnalysisTests/analysis/test_align.py +++ b/testsuite/MDAnalysisTests/analysis/test_align.py @@ -838,8 +838,10 @@ def test_alignto_rotates_velocities_and_forces(): ) reference.atoms.positions = ref_pos - angle = np.pi / 3 - known_R = transformations.rotation_matrix(angle, [0, 0, 1])[:3, :3] + angle = 23 + known_R = transformations.rotation_matrix(np.deg2rad(angle), [-1, 2, -3])[ + :3, :3 + ] mobile.atoms.positions = np.dot(ref_pos, known_R.T) rng = np.random.RandomState(0) diff --git a/testsuite/MDAnalysisTests/core/test_atomgroup.py b/testsuite/MDAnalysisTests/core/test_atomgroup.py index 79e4efcfe45..c05406f1d3a 100644 --- a/testsuite/MDAnalysisTests/core/test_atomgroup.py +++ b/testsuite/MDAnalysisTests/core/test_atomgroup.py @@ -439,9 +439,11 @@ def test_rotate_velocities_forces(self): def test_rotate_no_velocities_forces_does_not_raise(self): u = mda.Universe.empty(2, trajectory=True) - u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) - R = transformations.rotation_matrix(1, [0, 0, 1])[:3, :3] + orig_pos = np.array([[1, 0, 0], [-1, 0, 0]]) + u.atoms.positions = orig_pos.copy() + R = transformations.rotation_matrix(1, [-1, 2, -3])[:3, :3] u.atoms.rotate(R) + assert_almost_equal(u.atoms.positions, np.dot(orig_pos, R.T)) def test_transform_rotation_only(self, u, coords): R = np.eye(3) diff --git a/testsuite/MDAnalysisTests/transformations/test_rotate.py b/testsuite/MDAnalysisTests/transformations/test_rotate.py index c003d6535b1..9b4e0824944 100644 --- a/testsuite/MDAnalysisTests/transformations/test_rotate.py +++ b/testsuite/MDAnalysisTests/transformations/test_rotate.py @@ -192,9 +192,9 @@ def test_rotateby_velocities_forces(): orig_v = ts.velocities.copy() orig_f = ts.forces.copy() - axis = [0, 0, 1] + axis = [-1, 2, -3] point = [0, 0, 0] - angle = 90 + angle = 23 matrix = rotation_matrix(np.deg2rad(angle), axis, point) rotation = matrix[:3, :3].T @@ -210,9 +210,18 @@ def test_rotateby_velocities_forces(): def test_rotateby_no_velocities_forces_does_not_raise(): u = mda.Universe.empty(2, trajectory=True) - u.atoms.positions = np.array([[1, 0, 0], [-1, 0, 0]]) + orig_pos = np.array([[1, 0, 0], [-1, 0, 0]]) + u.atoms.positions = orig_pos.copy() ts = u.trajectory.ts - rotateby(90, [0, 0, 1], point=[0, 0, 0])(ts) + angle = 5 + matrix = rotation_matrix( + np.deg2rad(angle), [-1, 2, -3], [0, 0, 0] + ) + rotation = matrix[:3, :3].T + transformed_ts = rotateby(angle, [-1, 2, -3], point=[0, 0, 0])(ts) + assert_array_almost_equal( + transformed_ts.positions, np.dot(orig_pos, rotation), decimal=6 + ) @pytest.mark.parametrize( From fd257eb46cea8b7e57f208fd4a3840c69ce041d5 Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Sun, 9 Aug 2026 22:41:19 +0530 Subject: [PATCH 4/7] reformat using black --- testsuite/MDAnalysisTests/transformations/test_rotate.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/testsuite/MDAnalysisTests/transformations/test_rotate.py b/testsuite/MDAnalysisTests/transformations/test_rotate.py index 9b4e0824944..2ff376a0194 100644 --- a/testsuite/MDAnalysisTests/transformations/test_rotate.py +++ b/testsuite/MDAnalysisTests/transformations/test_rotate.py @@ -214,9 +214,7 @@ def test_rotateby_no_velocities_forces_does_not_raise(): u.atoms.positions = orig_pos.copy() ts = u.trajectory.ts angle = 5 - matrix = rotation_matrix( - np.deg2rad(angle), [-1, 2, -3], [0, 0, 0] - ) + matrix = rotation_matrix(np.deg2rad(angle), [-1, 2, -3], [0, 0, 0]) rotation = matrix[:3, :3].T transformed_ts = rotateby(angle, [-1, 2, -3], point=[0, 0, 0])(ts) assert_array_almost_equal( From e99d6d95ed55c8a4caaba28cec7351f460403ac8 Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Mon, 10 Aug 2026 22:08:03 +0530 Subject: [PATCH 5/7] remove redundant tests --- testsuite/MDAnalysisTests/core/test_atomgroup.py | 8 -------- .../MDAnalysisTests/transformations/test_rotate.py | 14 -------------- 2 files changed, 22 deletions(-) diff --git a/testsuite/MDAnalysisTests/core/test_atomgroup.py b/testsuite/MDAnalysisTests/core/test_atomgroup.py index c05406f1d3a..08085113431 100644 --- a/testsuite/MDAnalysisTests/core/test_atomgroup.py +++ b/testsuite/MDAnalysisTests/core/test_atomgroup.py @@ -437,14 +437,6 @@ def test_rotate_velocities_forces(self): assert_almost_equal(u.atoms.velocities, np.dot(orig_v, R.T)) assert_almost_equal(u.atoms.forces, np.dot(orig_f, R.T)) - def test_rotate_no_velocities_forces_does_not_raise(self): - u = mda.Universe.empty(2, trajectory=True) - orig_pos = np.array([[1, 0, 0], [-1, 0, 0]]) - u.atoms.positions = orig_pos.copy() - R = transformations.rotation_matrix(1, [-1, 2, -3])[:3, :3] - u.atoms.rotate(R) - assert_almost_equal(u.atoms.positions, np.dot(orig_pos, R.T)) - def test_transform_rotation_only(self, u, coords): R = np.eye(3) u.atoms.rotate(R) diff --git a/testsuite/MDAnalysisTests/transformations/test_rotate.py b/testsuite/MDAnalysisTests/transformations/test_rotate.py index 2ff376a0194..1b55aa2e2b2 100644 --- a/testsuite/MDAnalysisTests/transformations/test_rotate.py +++ b/testsuite/MDAnalysisTests/transformations/test_rotate.py @@ -208,20 +208,6 @@ def test_rotateby_velocities_forces(): ) -def test_rotateby_no_velocities_forces_does_not_raise(): - u = mda.Universe.empty(2, trajectory=True) - orig_pos = np.array([[1, 0, 0], [-1, 0, 0]]) - u.atoms.positions = orig_pos.copy() - ts = u.trajectory.ts - angle = 5 - matrix = rotation_matrix(np.deg2rad(angle), [-1, 2, -3], [0, 0, 0]) - rotation = matrix[:3, :3].T - transformed_ts = rotateby(angle, [-1, 2, -3], point=[0, 0, 0])(ts) - assert_array_almost_equal( - transformed_ts.positions, np.dot(orig_pos, rotation), decimal=6 - ) - - @pytest.mark.parametrize( "ag", ( From 5cd7c4e7ca98ebdc8076d3aaf9298fd76cad6764 Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Mon, 10 Aug 2026 22:41:15 +0530 Subject: [PATCH 6/7] update docs --- package/MDAnalysis/core/groups.py | 8 ++++++-- package/MDAnalysis/transformations/rotate.py | 2 ++ 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/package/MDAnalysis/core/groups.py b/package/MDAnalysis/core/groups.py index c30ede06833..2bd966064c5 100644 --- a/package/MDAnalysis/core/groups.py +++ b/package/MDAnalysis/core/groups.py @@ -1601,8 +1601,9 @@ def rotate(self, R, point=(0, 0, 0)): rotateby : rotate around given axis and angle MDAnalysis.lib.transformations : module of all coordinate transforms - .. versionchanged: 2.11.0 - Also rotate velocities and forces if present in Timestep. + + .. versionchanged:: 2.11.0 + Also rotate velocities and forces if present in Timestep. """ R = np.asarray(R) point = np.asarray(point) @@ -1662,6 +1663,9 @@ def rotateby(self, angle, axis, point=None): MDAnalysis.lib.transformations.rotation_matrix : calculate :math:`\mathsf{R}` + + .. versionchanged:: 2.11.0 + Also rotate velocities and forces if present in Timestep. """ alpha = np.radians(angle) axis = np.asarray(axis) diff --git a/package/MDAnalysis/transformations/rotate.py b/package/MDAnalysis/transformations/rotate.py index c23152cd790..23ee490d08a 100644 --- a/package/MDAnalysis/transformations/rotate.py +++ b/package/MDAnalysis/transformations/rotate.py @@ -121,6 +121,8 @@ class rotateby(TransformationBase): .. versionchanged:: 2.0.0 The transformation was changed to inherit from the base class for limiting threads and checking if it can be used in parallel analysis. + .. versionchanged:: 2.11.0 + Also rotate velocities and forces if present in Timestep. """ def __init__( From daa0a7cfeb92a946b2523942b6b7f478f7cd4ac2 Mon Sep 17 00:00:00 2001 From: ParthUppal523 Date: Mon, 10 Aug 2026 22:53:21 +0530 Subject: [PATCH 7/7] reformat using black --- package/MDAnalysis/core/groups.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/package/MDAnalysis/core/groups.py b/package/MDAnalysis/core/groups.py index 2bd966064c5..af084926cd4 100644 --- a/package/MDAnalysis/core/groups.py +++ b/package/MDAnalysis/core/groups.py @@ -1601,7 +1601,7 @@ def rotate(self, R, point=(0, 0, 0)): rotateby : rotate around given axis and angle MDAnalysis.lib.transformations : module of all coordinate transforms - + .. versionchanged:: 2.11.0 Also rotate velocities and forces if present in Timestep. """ @@ -1663,7 +1663,7 @@ def rotateby(self, angle, axis, point=None): MDAnalysis.lib.transformations.rotation_matrix : calculate :math:`\mathsf{R}` - + .. versionchanged:: 2.11.0 Also rotate velocities and forces if present in Timestep. """