Add functionality for parsing CMAP energy records

src/BioSimSpace/Process/_amber.py

@@ -1413,6 +1413,74 @@ def getCurrentDihedralEnergy(self, time_series=False, region=0, soft_core=False)
             time_series=time_series, region=region, soft_core=soft_core, block=False
         )
 
+    def getCMAPEnergy(self, time_series=False, region=0, soft_core=False, block="AUTO"):
+        """
+        Get the CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        region : int
+            The region to which the record corresponds. There will only be more
+            than one region for FreeEnergy protocols, where 1 indicates the second
+            TI region.
+
+        soft_core : bool
+            Whether to get the record for the soft-core part of the system for the
+            chosen region.
+
+        block : bool
+            Whether to block until the process has finished running.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+           The CMAP energy.
+        """
+        from .. import Units as _Units
+
+        return self.getRecord(
+            "CMAP",
+            time_series=time_series,
+            unit=_Units.Energy.kcal_per_mol,
+            region=region,
+            soft_core=soft_core,
+            block=block,
+        )
+
+    def getCurrentCMAPEnergy(self, time_series=False, region=0, soft_core=False):
+        """
+        Get the current CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        region : int
+            The region to which the record corresponds. There will only be more
+            than one region for FreeEnergy protocols, where 1 indicates the second
+            TI region.
+
+        soft_core : bool
+            Whether to get the record for the soft-core part of the system for the
+            chosen region.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+           The CMAP energy.
+        """
+        return self.getCMAPEnergy(
+            time_series=time_series, region=region, soft_core=soft_core, block=False
+        )
+
     def getElectrostaticEnergy(
         self, time_series=False, region=0, soft_core=False, block="AUTO"
     ):

src/BioSimSpace/Process/_gromacs.py

@@ -1382,6 +1382,47 @@ def getCurrentImproperEnergy(self, time_series=False):
         """
         return self.getImproperEnergy(time_series, block=False)
 
+    def getCMAPEnergy(self, time_series=False, block="AUTO"):
+        """
+        Get the CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        block : bool
+            Whether to block until the process has finished running.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+            The CMAP energy.
+        """
+        from .. import Units as _Units
+
+        return self.getRecord("CMAPDIH", time_series, _Units.Energy.kj_per_mol, block)
+
+    def getCurrentCMAPEnergy(self, time_series=False):
+        """
+        Get the current CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+            The CMAP energy.
+        """
+        return self.getCMAPEnergy(time_series, block=False)
+
     def getLennardJones14(self, time_series=False, block="AUTO"):
         """
         Get the Lennard-Jones energy between atoms 1 and 4.

src/BioSimSpace/Sandpit/Exscientia/Process/_amber.py

@@ -1497,6 +1497,74 @@ def getCurrentDihedralEnergy(self, time_series=False, region=0, soft_core=False)
             time_series=time_series, region=region, soft_core=soft_core, block=False
         )
 
+    def getCMAPEnergy(self, time_series=False, region=0, soft_core=False, block="AUTO"):
+        """
+        Get the CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        region : int
+            The region to which the record corresponds. There will only be more
+            than one region for FreeEnergy protocols, where 1 indicates the second
+            TI region.
+
+        soft_core : bool
+            Whether to get the record for the soft-core part of the system for the
+            chosen region.
+
+        block : bool
+            Whether to block until the process has finished running.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+           The CMAP energy.
+        """
+        from .. import Units as _Units
+
+        return self.getRecord(
+            "CMAP",
+            time_series=time_series,
+            unit=_Units.Energy.kcal_per_mol,
+            region=region,
+            soft_core=soft_core,
+            block=block,
+        )
+
+    def getCurrentCMAPEnergy(self, time_series=False, region=0, soft_core=False):
+        """
+        Get the current CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        region : int
+            The region to which the record corresponds. There will only be more
+            than one region for FreeEnergy protocols, where 1 indicates the second
+            TI region.
+
+        soft_core : bool
+            Whether to get the record for the soft-core part of the system for the
+            chosen region.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+           The CMAP energy.
+        """
+        return self.getCMAPEnergy(
+            time_series=time_series, region=region, soft_core=soft_core, block=False
+        )
+
     def getElectrostaticEnergy(
         self, time_series=False, region=0, soft_core=False, block="AUTO"
     ):

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

@@ -1440,6 +1440,47 @@ def getCurrentImproperEnergy(self, time_series=False):
         """
         return self.getImproperEnergy(time_series, block=False)
 
+    def getCMAPEnergy(self, time_series=False, block="AUTO"):
+        """
+        Get the CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        block : bool
+            Whether to block until the process has finished running.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+            The CMAP energy.
+        """
+        from .. import Units as _Units
+
+        return self.getRecord("CMAPDIH", time_series, _Units.Energy.kj_per_mol, block)
+
+    def getCurrentCMAPEnergy(self, time_series=False):
+        """
+        Get the current CMAP energy.
+
+        Parameters
+        ----------
+
+        time_series : bool
+            Whether to return a list of time series records.
+
+        Returns
+        -------
+
+        energy : :class:`Energy <BioSimSpace.Types.Energy>`
+            The CMAP energy.
+        """
+        return self.getCMAPEnergy(time_series, block=False)
+
     def getLennardJones14(self, time_series=False, block="AUTO"):
         """
         Get the Lennard-Jones energy between atoms 1 and 4.

tests/Process/test_single_point_energy.py

@@ -12,6 +12,16 @@ def ubiquitin_system():
     )
 
 
+@pytest.fixture(scope="module")
+def cmap_system():
+    """An ff19SB system, which carries CMAP backbone correction terms."""
+    import sire as sr
+
+    return BSS._SireWrappers.System(
+        sr.load_test_files("zero_k_torsions.gro", "zero_k_torsions.top")._system
+    )
+
+
 @pytest.mark.skipif(
     has_amber is False or has_gromacs is False,
     reason="Requires that both AMBER and GROMACS are installed.",
@@ -100,3 +110,42 @@ def test_amber_gromacs_triclinic(ubiquitin_system):
     nrg_amb = process_amb.getDihedralEnergy().kj_per_mol().value()
     nrg_gmx = process_gmx.getDihedralEnergy().kj_per_mol().value()
     assert nrg_amb == pytest.approx(nrg_gmx, rel=1e-2)
+
+
+@pytest.mark.skipif(
+    has_amber is False or has_gromacs is False,
+    reason="Requires that both AMBER and GROMACS are installed.",
+)
+def test_amber_gromacs_cmap(cmap_system):
+    """Single point CMAP energy comparison between AMBER and GROMACS."""
+
+    # Create a single-step minimisation protocol.
+    protocol = BSS.Protocol.Minimisation(steps=1)
+
+    # Create a process to run with AMBER.
+    process_amb = BSS.Process.Amber(cmap_system, protocol)
+
+    # Create a process to run with GROMACS.
+    process_gmx = BSS.Process.Gromacs(
+        cmap_system, protocol, extra_options={"nsteps": 0}
+    )
+
+    # Run the AMBER process and wait for it to finish.
+    process_amb.start()
+    process_amb.wait()
+
+    # Run the GROMACS process and wait for it to finish.
+    process_gmx.start()
+    process_gmx.wait()
+
+    # Compare CMAP energies. (In kJ / mol)
+    nrg_amb = process_amb.getCMAPEnergy()
+    nrg_gmx = process_gmx.getCMAPEnergy()
+
+    # The comparison is meaningless if either engine didn't report the term.
+    assert nrg_amb is not None
+    assert nrg_gmx is not None
+
+    assert nrg_amb.kj_per_mol().value() == pytest.approx(
+        nrg_gmx.kj_per_mol().value(), rel=1e-2
+    )

tests/Sandpit/Exscientia/Process/test_single_point_energy.py

@@ -12,6 +12,16 @@ def system():
     )
 
 
+@pytest.fixture(scope="session")
+def cmap_system():
+    """An ff19SB system, which carries CMAP backbone correction terms."""
+    import sire as sr
+
+    return BSS._SireWrappers.System(
+        sr.load_test_files("zero_k_torsions.gro", "zero_k_torsions.top")._system
+    )
+
+
 @pytest.mark.skipif(
     has_amber is False or has_gromacs is False or has_pyarrow is False,
     reason="Requires that AMBER, GROMACS, and pyarrow are installed.",
@@ -96,3 +106,42 @@ def test_amber_gromacs_triclinic(system):
     nrg_amb = process_amb.getDihedralEnergy().kj_per_mol().value()
     nrg_gmx = process_gmx.getDihedralEnergy().kj_per_mol().value()
     assert nrg_amb == pytest.approx(nrg_gmx, rel=1e-2)
+
+
+@pytest.mark.skipif(
+    has_amber is False or has_gromacs is False or has_pyarrow is False,
+    reason="Requires that AMBER, GROMACS, and pyarrow are installed.",
+)
+def test_amber_gromacs_cmap(cmap_system):
+    """Single point CMAP energy comparison between AMBER and GROMACS."""
+
+    # Create a single-step minimisation protocol.
+    protocol = BSS.Protocol.Minimisation(steps=1)
+
+    # Create a process to run with AMBER.
+    process_amb = BSS.Process.Amber(cmap_system, protocol)
+
+    # Create a process to run with GROMACS.
+    process_gmx = BSS.Process.Gromacs(
+        cmap_system, protocol, extra_options={"nsteps": 0}
+    )
+
+    # Run the AMBER process and wait for it to finish.
+    process_amb.start()
+    process_amb.wait()
+
+    # Run the GROMACS process and wait for it to finish.
+    process_gmx.start()
+    process_gmx.wait()
+
+    # Compare CMAP energies. (In kJ / mol)
+    nrg_amb = process_amb.getCMAPEnergy()
+    nrg_gmx = process_gmx.getCMAPEnergy()
+
+    # The comparison is meaningless if either engine didn't report the term.
+    assert nrg_amb is not None
+    assert nrg_gmx is not None
+
+    assert nrg_amb.kj_per_mol().value() == pytest.approx(
+        nrg_gmx.kj_per_mol().value(), rel=1e-2
+    )