Merge pull request #1045 from openforcefield/fix-1033

Fix charge ordering in `.prmtop` files

docs/releasehistory.md

@@ -11,6 +11,10 @@ Dates are given in YYYY-MM-DD format.
 
 Please note that all releases prior to a version 1.0.0 are considered pre-releases and many API changes will come before a stable release.
 
+## 0.3.x
+
+* #1045 Fixes a bug in which charges were sometimes misordered in `.prmtop` files.
+
 ## 0.3.30 - 2024-08
 
 * #1039 Updates support of "cutoff" electrostatics in `.to_openmm` to better reflect what OpenMM supports. Set `"reaction-field"` to force the use of `CutoffPeriodic`, provided the vdW and electrostatic cutoff distances match. The potential/method `"cutoff"` is no longer supported but may be re-added in the future.

openff/interchange/_tests/test_issues.py

@@ -114,3 +114,23 @@ def test_issue_1031(monkeypatch):
     # check a few atom names to ensure these didn't end up being empty sets
     for atom_name in ("NE2", "H3", "HA", "CH3", "CA", "CB", "CE1"):
         assert atom_name in openff_atom_names
+
+
+def test_issue_1033(water, sage):
+    """Test issue/PR 1033, in which charges were not ordered topologically."""
+    pytest.importorskip("parmed")
+
+    interchange = sage.create_interchange(
+        Topology.from_molecules(3 * [water]),
+    )
+
+    interchange.to_prmtop("test.prmtop")
+
+    structure = parmed.load_file("test.prmtop")
+
+    # atoms should be ordered OHH OHH OHH, charges should be as well
+    assert 3 * [
+        (8, -0.834),
+        (1, 0.417),
+        (1, 0.417),
+    ] == [(atom.atomic_number, atom.charge) for atom in structure.atoms]

openff/interchange/interop/amber/export/_export.py

@@ -22,7 +22,7 @@
     UnsupportedExportError,
     UnsupportedMixingRuleError,
 )
-from openff.interchange.models import PotentialKey
+from openff.interchange.models import PotentialKey, TopologyKey
 
 
 def _write_text_blob(file, blob):
@@ -514,9 +514,15 @@ def to_prmtop(interchange: "Interchange", file_path: Path | str):
         _write_text_blob(prmtop, text_blob)
 
         prmtop.write("%FLAG CHARGE\n" "%FORMAT(5E16.8)\n")
+        electro = interchange["Electrostatics"]
         charges = [
-            charge.m_as(unit.e) * AMBER_COULOMBS_CONSTANT
-            for charge in interchange["Electrostatics"].charges.values()
+            charge * AMBER_COULOMBS_CONSTANT
+            for charge in [
+                electro.charges[TopologyKey(atom_indices=(index,))].m_as(
+                    unit.elementary_charge,
+                )
+                for index in range(interchange.topology.n_atoms)
+            ]
         ]
         text_blob = "".join([f"{val:16.8E}" for val in charges])
         _write_text_blob(prmtop, text_blob)