core: bond_breakage: handle ActionType::NONE as a no-op (bug-sweep #33)

src/core/bond_breakage/bond_breakage.cpp

@@ -79,7 +79,6 @@ static ActionSet actions_for_breakage(CellStructure const &cell_structure,
     return bond_partners[1];
   }; // optional for second partner engaged
 
-  // Handle different action types
   if (spec.action_type == ActionType::DELETE_BOND) {
     if (is_angle_bond(e.bond_partners)) {
       return {DeleteAngleBond{e.particle_id,
@@ -88,54 +87,59 @@ static ActionSet actions_for_breakage(CellStructure const &cell_structure,
     }
     return {DeleteBond{e.particle_id, *(e.bond_partners[0]), e.bond_type}};
   }
+
 #ifdef ESPRESSO_VIRTUAL_SITES_RELATIVE
   // revert bind at point of collision for pair bonds
-  if (spec.action_type == ActionType::REVERT_BIND_AT_POINT_OF_COLLISION and
-      not is_angle_bond(e.bond_partners)) {
-    // We need to find the base particles for the two virtual sites
-    // between which the bond broke.
-    auto p1 = cell_structure.get_local_particle(e.particle_id);
-    auto p2 = cell_structure.get_local_particle(*(e.bond_partners[0]));
-    if (p1 and p2) {
-      if (not p1->is_virtual() or not p2->is_virtual()) {
-        runtimeErrorMsg() << "The REVERT_BIND_AT_POINT_OF_COLLISION bond "
-                             "breakage action has to be configured for the "
-                             "bond on the virtual site. Encountered a particle "
-                             "that is not virtual.";
-        return {};
-      }
+  if (spec.action_type == ActionType::REVERT_BIND_AT_POINT_OF_COLLISION) {
+    if (not is_angle_bond(e.bond_partners)) {
+      // We need to find the base particles for the two virtual sites
+      // between which the bond broke.
+      auto p1 = cell_structure.get_local_particle(e.particle_id);
+      auto p2 = cell_structure.get_local_particle(*(e.bond_partners[0]));
+      if (p1 and p2) {
+        if (not p1->is_virtual() or not p2->is_virtual()) {
+          runtimeErrorMsg()
+              << "The REVERT_BIND_AT_POINT_OF_COLLISION bond "
+                 "breakage action has to be configured for the "
+                 "bond on the virtual site. Encountered a particle "
+                 "that is not virtual.";
+          return {};
+        }
 
-      return {
-          // Bond between virtual sites
-          DeleteBond{e.particle_id, *(e.bond_partners[0]), e.bond_type},
-          // Bond between base particles. We do not know, on which of these
-          // the bond is defined, since bonds are stored only on one partner
-          DeleteAllBonds{p1->vs_relative().to_particle_id,
-                         p2->vs_relative().to_particle_id},
-          DeleteAllBonds{p2->vs_relative().to_particle_id,
-                         p1->vs_relative().to_particle_id},
-      };
-    }
-  } else {
-    // revert bind at point of collision for angle bonds
-    auto vs = cell_structure.get_local_particle(e.particle_id);
-    auto p1 = cell_structure.get_local_particle(*(e.bond_partners[0]));
-    auto p2 = cell_structure.get_local_particle(*(e.bond_partners[1]));
-    if (p1 and p2) {
-      if (not vs->is_virtual()) {
-        runtimeErrorMsg() << "The REVERT_BIND_AT_POINT_OF_COLLISION bond "
-                             "breakage action has to be configured for the "
-                             "bond on the virtual site. Encountered a particle "
-                             "that is not virtual.";
-        return {};
+        return {
+            // Bond between virtual sites
+            DeleteBond{e.particle_id, *(e.bond_partners[0]), e.bond_type},
+            // Bond between base particles. We do not know, on which of these
+            // the bond is defined, since bonds are stored only on one partner
+            DeleteAllBonds{p1->vs_relative().to_particle_id,
+                           p2->vs_relative().to_particle_id},
+            DeleteAllBonds{p2->vs_relative().to_particle_id,
+                           p1->vs_relative().to_particle_id},
+        };
       }
+    } else {
+      // revert bind at point of collision for angle bonds
+      auto vs = cell_structure.get_local_particle(e.particle_id);
+      auto p1 = cell_structure.get_local_particle(*(e.bond_partners[0]));
+      auto p2 = cell_structure.get_local_particle(*(e.bond_partners[1]));
+      if (p1 and p2) {
+        if (not vs->is_virtual()) {
+          runtimeErrorMsg()
+              << "The REVERT_BIND_AT_POINT_OF_COLLISION bond "
+                 "breakage action has to be configured for the "
+                 "bond on the virtual site. Encountered a particle "
+                 "that is not virtual.";
+          return {};
+        }
 
-      return {// Angle bond on the virtual site
-              DeleteAngleBond{e.particle_id, {p1->id(), p2->id()}, e.bond_type},
-              // Bond between base particles. We do not know, on which of these
-              // the bond is defined, since bonds are stored only on one partner
-              DeleteAllBonds{p1->id(), p2->id()},
-              DeleteAllBonds{p2->id(), p1->id()}};
+        return {
+            // Angle bond on the virtual site
+            DeleteAngleBond{e.particle_id, {p1->id(), p2->id()}, e.bond_type},
+            // Bond between base particles. We do not know, on which of these
+            // the bond is defined, since bonds are stored only on one partner
+            DeleteAllBonds{p1->id(), p2->id()},
+            DeleteAllBonds{p2->id(), p1->id()}};
+      }
     }
   }
 #endif // ESPRESSO_VIRTUAL_SITES_RELATIVE

testsuite/python/bond_breakage.py

@@ -253,18 +253,15 @@ def count_bonds(self, pairs):
         return bonds_count
 
     def setUp(self):
-
-        box_vol = self.system.box_l[0]**3.
         phi = 0.4
-
         r = 1.
-        solid_vol = phi * box_vol
+        solid_vol = phi * self.system.volume()
         part_vol = 4 / 3 * np.pi * r**3
         part_num = int(solid_vol / part_vol)
 
         np.random.seed(seed=678)
         for i in range(part_num):
-            pos = np.random.rand(3) * self.system.box_l[0]
+            pos = np.random.rand(3) * self.system.box_l
             self.system.part.add(pos=pos)
 
         self.system.non_bonded_inter[0, 0].lennard_jones.set_params(
@@ -360,5 +357,42 @@ def test_vs_bonds(self):
         np.testing.assert_equal(bonds_dist, bonds_count)
 
 
+class BreakageAPI(BondBreakageCommon, ut.TestCase):
+
+    def setUp(self):
+        self.system.box_l = 3 * [20]
+        self.system.min_global_cut = 0.1
+        self.system.time_step = 0.01
+        self.system.cell_system.skin = 0.4
+
+    def tearDown(self):
+        self.system.part.clear()
+        self.system.bonded_inter.clear()
+
+    def test_bond_deletion(self):
+        p1 = self.system.part.add(pos=[0., 0., 0.], v=[0., 0., 0.])
+        p2 = self.system.part.add(pos=[0., 0., 1.], v=[0., 0., 1.])
+        bond = espressomd.interactions.HarmonicBond(k=1., r_0=1., r_cut=1.2)
+        self.system.bonded_inter.add(bond)
+        p1.bonds = ((bond, p2))
+        self.system.bond_breakage[bond] = BreakageSpec(
+            breakage_length=1.1, action_type="delete_bond")
+        self.system.integrator.run(100)
+        assert np.linalg.norm(p2.pos - p1.pos) > bond.r_cut
+        self.assertEqual(len(p1.bonds), 0)
+
+    def test_none(self):
+        p1 = self.system.part.add(pos=[0., 0., 0.], v=[0., 0., 0.])
+        p2 = self.system.part.add(pos=[0., 0., 1.], v=[0., 0., 1.])
+        bond = espressomd.interactions.HarmonicBond(k=1., r_0=1., r_cut=1.2)
+        self.system.bonded_inter.add(bond)
+        p1.bonds = ((bond, p2))
+        self.system.bond_breakage[bond] = BreakageSpec(
+            breakage_length=1.1, action_type="none")
+        self.system.integrator.run(100)
+        assert np.linalg.norm(p2.pos - p1.pos) > bond.r_cut
+        self.assertEqual(len(p1.bonds), 1)
+
+
 if __name__ == "__main__":
     ut.main()