Fixing vacuum translations in GROMACS (#57)

python/BioSimSpace/Sandpit/Exscientia/Process/_gromacs.py

@@ -47,6 +47,7 @@
 from sire.legacy import Units as _SireUnits
 from sire.legacy import Vol as _SireVol
 
+from ..Units.Length import angstrom
 from .._Utils import _assert_imported, _have_imported, _try_import
 
 # alchemlyb isn't available on all variants of Python that we support, so we
@@ -2704,8 +2705,22 @@ def _getFinalFrame(self):
                     )
 
             # If this is a vacuum simulation, then translate the centre of mass
-            # of the system back to the origin.
-            if not space_prop in old_system._sire_object.propertyKeys():
+            # of the system back to the middle of the box to preserve PBC.
+            if old_system.getBox() == (None, None):
+                try:
+                    old_box = old_system._sire_object.property(space_prop)
+                except:
+                    old_box = None
+                box = _SireVol.PeriodicBox(_SireMaths.Vector(9999, 9999, 9999))
+                old_system._sire_object.setProperty(space_prop, box)
+                com = [angstrom * 9999 / 2 for _ in range(3)]
+                old_system.translate([x for x in com])
+                old_system._sire_object.make_whole()
+                old_system.translate([-x for x in com])
+                if old_box is None:
+                    old_system._sire_object.removeProperty(space_prop)
+                else:
+                    old_system._sire_object.setProperty(space_prop, old_box)
                 com = old_system._getCenterOfMass()
                 old_system.translate([-x for x in com])
 
@@ -2821,9 +2836,23 @@ def _getFrame(self, time):
                             )
 
                 # If this is a vacuum simulation, then translate the centre of mass
-                # of the system back to the origin.
-                if not space_prop in old_system._sire_object.propertyKeys():
-                    com = new_system._getCenterOfMass()
+                # of the system back to the middle of the box to preserve PBC.
+                if old_system.getBox() == (None, None):
+                    try:
+                        old_box = old_system._sire_object.property(space_prop)
+                    except:
+                        old_box = None
+                    box = _SireVol.PeriodicBox(_SireMaths.Vector(9999, 9999, 9999))
+                    old_system._sire_object.setProperty(space_prop, box)
+                    com = [angstrom * 9999 / 2 for _ in range(3)]
+                    old_system.translate([x for x in com])
+                    old_system._sire_object.make_whole()
+                    old_system.translate([-x for x in com])
+                    if old_box is None:
+                        old_system._sire_object.removeProperty(space_prop)
+                    else:
+                        old_system._sire_object.setProperty(space_prop, old_box)
+                    com = old_system._getCenterOfMass()
                     old_system.translate([-x for x in com])
 
                 return old_system