fix: plumed_mod updating

Co-authored-by: Eric Hartmann <hartmaec@rh05659.villa-bosch.de>

plugins/src/homolysis/reaction.py

@@ -59,12 +59,9 @@ def get_recipe_collection(self, files: TaskFiles):
             # averaging distances works here because we typically have
             # one conformational state per calculation
 
-            # converto to zero-base
-            atomids = [i - 1 for i in list(atomids)]
-
             recipes.append(
                 Recipe(
-                    recipe_steps=[Break(*list(atomids))],
+                    recipe_steps=[Break(atom_id_1 = atomids[0], atom_id_2 = atomids[1])],
                     rates=[*k_avg],
                     timespans=[(distances["time"][0], distances["time"][-1])],
                 )

src/kimmdy/parsing.py

@@ -385,7 +385,7 @@ def read_plumed(path: Path) -> dict:
                     {
                         "id": d[0].strip(":"),
                         "keyword": d[1],
-                        "atoms": [int(x) for x in d[2].strip("ATOMS=").split(",")],
+                        "atoms": [str(x) for x in d[2].strip("ATOMS=").split(",")],
                     }
                 )
             elif "PRINT" in l[:5]:

src/kimmdy/utils.py

@@ -39,20 +39,19 @@ def increment_logfile(f: Path) -> Path:
 
 def get_atominfo_from_plumedid(
     plumedid: str, plumed: dict, top: dict
-) -> tuple[frozenset, frozenset]:
+) -> tuple[frozenset[str], list[str]]:
     """returns atomtypes for a plumedid with information from the plumed and topology file"""
     lookup_atomid_plumedid = {
         entry["id"]: frozenset(entry["atoms"]) for entry in plumed["distances"]
     }
     atoms = get_protein_section(top, "atoms")
     if not atoms:
         raise ValueError("Could not find atoms in topology file")
-    lookup_atomtype_atomid = {int(atom[0]): atom[1] for atom in atoms}
-    atomids = lookup_atomid_plumedid[plumedid]
-    atomids_list = list(atomids)
+    lookup_atomtype_atomid = {str(atom[0]): atom[1] for atom in atoms}
+    atomids = sorted(lookup_atomid_plumedid[plumedid],key=int)
     atomtypes_list = [
-        lookup_atomtype_atomid[atomids_list[0]],
-        lookup_atomtype_atomid[atomids_list[1]],
+        lookup_atomtype_atomid[atomids[0]],
+        lookup_atomtype_atomid[atomids[1]],
     ]
     atomtypes = frozenset(atomtypes_list)
     logging.debug(f"Found atomtypes {atomtypes} for plumedid {plumedid}.")

tests/test_changemanager.py

@@ -4,7 +4,7 @@
 
 from kimmdy.reaction import Break, Bind, Move, RecipeStep
 from kimmdy.parsing import read_plumed, read_top
-from kimmdy.changemanager import break_bond_plumed, modify_coords
+from kimmdy.changemanager import break_bond_plumed, modify_plumed,modify_coords
 from kimmdy.topology.topology import Topology
 from conftest import SlimFiles
 import os
@@ -32,12 +32,25 @@ def test_plumed_break(tmpdir):
     plumed = read_plumed(tmpdir / "plumed_nat.dat")
     plumed_break_ref = read_plumed(tmpdir / "plumed_break29-35.dat")
 
-    breakpair = (29, 35)
+    breakpair = ('29', '35')
     plumed_break = break_bond_plumed(plumed, breakpair, Path("distances.dat"))
 
     assert plumed_break["distances"] == plumed_break_ref["distances"]
     assert plumed_break["prints"] == plumed_break_ref["prints"]
 
+def test_plumed_modify(tmpdir):
+    plumeddat:Path = tmpdir / "plumed_nat.dat"
+    newplumeddat:Path = tmpdir / "plumed_test.dat"
+    recipe_steps = [Break(28,34)]
+    plumeddist: Path = Path("distances.dat")
+
+    modify_plumed(recipe_steps,plumeddat,newplumeddat,plumeddist)
+
+    plumed_break_ref = read_plumed(tmpdir / "plumed_break29-35.dat")
+    plumed_break_test = read_plumed(newplumeddat)
+
+    assert plumed_break_test["distances"] == plumed_break_ref["distances"]
+    assert plumed_break_test["prints"] == plumed_break_ref["prints"]
 
 def test_modify_coords_break(tmpdir):
     steps: list[RecipeStep] = [Break(28, 34)]

tests/test_homolysis.py

@@ -48,7 +48,7 @@ def test_lookup_atominfo(homolysis_files):
         "d1", homolysis_files["plumed"], homolysis_files["top"]
     )
     assert atomtypes == frozenset(["N", "CT"])
-    assert atomid == frozenset([7, 9])
+    assert atomid == ['7', '9']
 
 
 def test_lookup_bondprm(homolysis_files):