From deb3e3172068f18ff1ddedafdbc60a10b8a7f757 Mon Sep 17 00:00:00 2001 From: Rudolf Weeber Date: Mon, 15 Jun 2026 11:46:17 +0200 Subject: [PATCH] analysis: accept p_type=-1 for center_of_mass and angular_momentum (bug-sweep #50) The script-interface handlers for center_of_mass and angular_momentum unconditionally called check_particle_type(p_type), which throws for any p_type < 0. This rejected the documented and core-supported p_type == -1 sentinel ("all non-virtual particles"), so system.analysis.angular_momentum (p_type=-1) and center_of_mass(p_type=-1) raised "Particle type -1 does not exist" instead of computing the sum over all particles. The core functions center_of_mass and angular_momentum already implement the -1 case ((p.type() == p_type or p_type == -1)). Guard the validation so the documented sentinel bypasses check_particle_type while explicit types are still validated. The scope is intentionally limited to these two methods: min_dist and distribution still (correctly) reject -1. Also align the center_of_mass docstring with angular_momentum's already-correct text mentioning the -1 sentinel, and add regression tests asserting the -1 result equals the sum over all non-virtual particles (distinct from the type-0-only result, so a wrong-subset regression is also caught). Co-Authored-By: Claude Opus 4.8 --- src/python/espressomd/analyze.py | 3 ++- src/script_interface/analysis/Analysis.cpp | 14 +++++++++++-- testsuite/python/analyze_mass_related.py | 24 ++++++++++++++++++++++ 3 files changed, 38 insertions(+), 3 deletions(-) diff --git a/src/python/espressomd/analyze.py b/src/python/espressomd/analyze.py index 693267c643..ca2d4c106c 100644 --- a/src/python/espressomd/analyze.py +++ b/src/python/espressomd/analyze.py @@ -111,7 +111,8 @@ class Analysis(ScriptInterfaceHelper): ---------- p_type : :obj:`int` Particle :attr:`~espressomd.particle_data.ParticleHandle.type` - for which to calculate the center of mass. + for which to calculate the center of mass, or ``-1`` for all + particles. Returns ------- diff --git a/src/script_interface/analysis/Analysis.cpp b/src/script_interface/analysis/Analysis.cpp index d6aa626084..5398690585 100644 --- a/src/script_interface/analysis/Analysis.cpp +++ b/src/script_interface/analysis/Analysis.cpp @@ -167,13 +167,23 @@ Variant Analysis::do_call_method(std::string const &name, } if (name == "center_of_mass") { auto const p_type = get_value(parameters, "p_type"); - context()->parallel_try_catch([&]() { check_particle_type(p_type); }); + // p_type == -1 is the documented sentinel for "all (non-virtual) + // particles", which the core function handles directly; only validate + // explicit (non-sentinel) particle types. + if (p_type != -1) { + context()->parallel_try_catch([&]() { check_particle_type(p_type); }); + } auto const local = center_of_mass(get_system(), p_type); return mpi_reduce_sum(context()->get_comm(), local).as_vector(); } if (name == "angular_momentum") { auto const p_type = get_value(parameters, "p_type"); - context()->parallel_try_catch([&]() { check_particle_type(p_type); }); + // p_type == -1 is the documented sentinel for "all (non-virtual) + // particles", which the core function handles directly; only validate + // explicit (non-sentinel) particle types. + if (p_type != -1) { + context()->parallel_try_catch([&]() { check_particle_type(p_type); }); + } auto const local = angular_momentum(get_system(), p_type); return mpi_reduce_sum(context()->get_comm(), local).as_vector(); } diff --git a/testsuite/python/analyze_mass_related.py b/testsuite/python/analyze_mass_related.py index 66c8714eb9..b66c6b0575 100644 --- a/testsuite/python/analyze_mass_related.py +++ b/testsuite/python/analyze_mass_related.py @@ -97,6 +97,18 @@ def test_center_of_mass(self): np.testing.assert_allclose(com, com_ref) + def test_center_of_mass_all_particles(self): + # p_type=-1 is the documented sentinel for "all (non-virtual) + # particles", and must be distinct from the type-0-only result. + no_virtual = self.system.part.select(lambda p: not p.is_virtual()) + com_ref = np.zeros(3) + for p in no_virtual: + com_ref += p.pos * p.mass + com_ref /= np.sum(no_virtual.mass) + com = self.system.analysis.center_of_mass(p_type=-1) + + np.testing.assert_allclose(com, com_ref) + def test_galilei_transform(self): no_virtual = self.system.part.select(lambda p: not p.is_virtual()) @@ -120,6 +132,18 @@ def test_angularmomentum(self): am, self.system.analysis.angular_momentum(p_type=0)) + def test_angularmomentum_all_particles(self): + # p_type=-1 is the documented sentinel for "all (non-virtual) + # particles", and must be distinct from the type-0-only result. + no_virtual = self.system.part.select(lambda p: not p.is_virtual()) + am_ref = np.zeros(3) + for p in no_virtual: + am_ref += p.mass * np.cross(p.pos, p.v) + + np.testing.assert_allclose( + am_ref, + self.system.analysis.angular_momentum(p_type=-1)) + def test_kinetic_energy(self): no_virtual = self.system.part.select(lambda p: not p.is_virtual()) E_kin = 0.5 * \