[WIP] Handle NAGLCharges

devtools/conda-envs/docs_env.yaml

@@ -27,3 +27,4 @@ dependencies:
   - sphinx-notfound-page
   - pip:
       - git+https://github.com/openforcefield/openff-sphinx-theme.git@main
+      - git+https://github.com/openforcefield/openff-toolkit.git@naglcharges-handler

devtools/conda-envs/examples_env.yaml

@@ -38,3 +38,5 @@ dependencies:
   - pdbfixer
   - openeye-toolkits =2024.1.0
   - rich
+  - pip:
+      - git+https://github.com/openforcefield/openff-toolkit.git@naglcharges-handler

devtools/conda-envs/test_env.yaml

@@ -35,3 +35,6 @@ dependencies:
   - typing-extensions
   - types-setuptools
   - pandas-stubs
+  # Temporary testing dep
+  - pip:
+      - git+https://github.com/openforcefield/openff-toolkit.git@naglcharges-handler

devtools/conda-envs/test_not_py313.yaml

@@ -45,3 +45,5 @@ dependencies:
   - typing-extensions
   - types-setuptools
   - pandas-stubs
+  - pip:
+      - git+https://github.com/openforcefield/openff-toolkit.git@naglcharges-handler

openff/interchange/_tests/unit_tests/smirnoff/test_nonbonded.py

@@ -135,6 +135,33 @@ def test_toolkit_am1bcc_uses_elf10_if_oe_is_available(self, sage, hexane_diol):
             assert not uses_elf10
             numpy.testing.assert_allclose(partial_charges, assigned_charges)
 
+    def test_nagl_charge_assignment_matches_reference(self, hexane_diol):
+        from openff.toolkit.typing.engines.smirnoff import ForceField
+
+        from openff.interchange import Interchange
+
+        hexane_diol.assign_partial_charges("openff-gnn-am1bcc-0.1.0-rc.3.pt")
+        # Leave the ToolkitAM1BCC tag in openff-2.1.0 to ensure that the NAGLCharges handler takes precedence
+        ff = ForceField("openff-2.1.0.offxml")
+        ff.get_parameter_handler(
+            "NAGLCharges",
+            {
+                "model_file": "openff-gnn-am1bcc-0.1.0-rc.3.pt",
+                "version": "0.3",
+            },
+        )
+
+        interchange = Interchange.from_smirnoff(force_field=ff, topology=hexane_diol.to_topology())
+
+        assigned_charges_unitless = [v.m for v in interchange["Electrostatics"]._get_charges().values()]
+
+        expected_charges = hexane_diol.partial_charges
+        assert expected_charges is not None
+        assert expected_charges.units == unit.elementary_charge
+        assert not all(charge == 0 for charge in expected_charges.magnitude)
+        expected_charges_unitless = [v.m for v in expected_charges]
+        numpy.testing.assert_allclose(expected_charges_unitless, assigned_charges_unitless)
+
     @pytest.mark.skip(
         reason="Turn on if toolkit ever allows non-standard scale12/13/15",
     )

openff/interchange/smirnoff/_create.py

@@ -69,7 +69,7 @@ def _check_supported_handlers(force_field: ForceField):
     unsupported = list()
 
     for handler_name in force_field.registered_parameter_handlers:
-        if handler_name in {"ToolkitAM1BCC"}:
+        if handler_name in {"ToolkitAM1BCC", "NAGLCharges"}:
             continue
         if handler_name not in _SUPPORTED_PARAMETER_HANDLERS:
             unsupported.append(handler_name)
@@ -348,6 +348,7 @@ def _electrostatics(
                         force_field._parameter_handlers.get(name, None)
                         for name in [
                             "Electrostatics",
+                            "NAGLCharges",
                             "ChargeIncrementModel",
                             "ToolkitAM1BCC",
                             "LibraryCharges",

openff/interchange/smirnoff/_nonbonded.py

@@ -11,10 +11,12 @@
     ChargeIncrementModelHandler,
     ElectrostaticsHandler,
     LibraryChargeHandler,
+    NAGLChargesHandler,
     ParameterHandler,
     ToolkitAM1BCCHandler,
     vdWHandler,
 )
+from openff.toolkit.utils.exceptions import MissingPackageError
 from pydantic import Field, PrivateAttr, computed_field
 from typing_extensions import Self
 
@@ -48,6 +50,7 @@
 
 ElectrostaticsHandlerType = Union[
     ElectrostaticsHandler,
+    NAGLChargesHandler,
     ToolkitAM1BCCHandler,
     ChargeIncrementModelHandler,
     LibraryChargeHandler,
@@ -248,7 +251,7 @@
     * global settings for the electrostatic interactions such as the cutoff distance
       and the intramolecular scale factors.
     * partial charges which have been assigned by a ``ToolkitAM1BCC``,
-      ``LibraryCharges``, or a ``ChargeIncrementModel`` parameter
+      ``LibraryCharges``, ``NAGLCharges``, or a ``ChargeIncrementModel`` parameter
       handler.
     * charge corrections applied by a ``ChargeIncrementHandler``
 
@@ -284,6 +287,7 @@
         """Return a list of allowed types of ParameterHandler classes."""
         return [
             LibraryChargeHandler,
+            NAGLChargesHandler,
             ChargeIncrementModelHandler,
             ToolkitAM1BCCHandler,
             ElectrostaticsHandler,
@@ -364,6 +368,7 @@
 
                     if potential_key.associated_handler in (
                         "LibraryCharges",
+                        "NAGLChargesHandler",
                         "ToolkitAM1BCCHandler",
                         "molecules_with_preset_charges",
                         "ExternalSource",
@@ -430,7 +435,7 @@
         """
         Return the order in which parameter handlers take precedence when computing charges.
         """
-        return ["LibraryCharges", "ChargeIncrementModel", "ToolkitAM1BCC"]
+        return ["LibraryCharges", "NAGLCharges", "ChargeIncrementModel", "ToolkitAM1BCC"]
 
     @classmethod
     def create(
@@ -642,7 +647,7 @@
     @classmethod
     def _find_charge_model_matches(
         cls,
-        parameter_handler: ToolkitAM1BCCHandler | ChargeIncrementModelHandler,
+        parameter_handler: ToolkitAM1BCCHandler | ChargeIncrementModelHandler | NAGLChargesHandler,
         unique_molecule: Molecule,
     ) -> tuple[
         str,
@@ -659,8 +664,19 @@
 
         handler_name = parameter_handler.__class__.__name__
 
-        if handler_name == "ChargeIncrementModelHandler":
+        if handler_name == "NAGLChargesHandler":
+            from openff.toolkit.utils.toolkits import GLOBAL_TOOLKIT_REGISTRY
+
+            partial_charge_method = parameter_handler.model_file
+            if "NAGL" not in GLOBAL_TOOLKIT_REGISTRY.__repr__():
+                raise MissingPackageError(
+                    "The force field has a NAGLCharges section, but the NAGL software isn't "
+                    "present in GLOBAL_TOOLKIT_REGISTRY",
+                )
+
+        elif handler_name == "ChargeIncrementModelHandler":
             partial_charge_method = parameter_handler.partial_charge_method
+
         elif handler_name == "ToolkitAM1BCCHandler":
             from openff.toolkit.utils.toolkits import GLOBAL_TOOLKIT_REGISTRY
 
@@ -673,9 +689,12 @@
         else:
             raise InvalidParameterHandlerError(
                 f"Encountered unknown handler of type {type(parameter_handler)} where only "
-                "ToolkitAM1BCCHandler or ChargeIncrementModelHandler are expected.",
+                "ToolkitAM1BCCHandler, NAGLChargesHandler, or ChargeIncrementModelHandler are expected.",
             )
 
+        # Comment for reviewer: Could get fancy here and force ToolkitAM1BCCHandler calls to go to
+        # AmberTools/OpenEyeToolkitWrapper, and NAGLChargesHandler to go to NAGLToolkitWrapper, but my initial
+        # judgement here is that this isn't worth the complexity. Happy to change this if that's the case.
         partial_charges = cls._compute_partial_charges(
             unique_molecule,
             unique_molecule.to_smiles(
@@ -741,7 +760,7 @@
                     unique_molecule,
                 )
 
-            if handler_type in ["ToolkitAM1BCC", "ChargeIncrementModel"]:
+            if handler_type in ["ToolkitAM1BCC", "ChargeIncrementModel", "NAGLCharges"]:
                 (
                     partial_charge_method,
                     am1_matches,
@@ -922,6 +941,12 @@
                                         f"{new_key.extras['partial_charge_method']}, "
                                         f"applied to topology atom index {topology_atom_index}",
                                     )
+                                elif new_key.extras["handler"] == "NAGLChargesHandler":
+                                    logger.info(
+                                        "Charge section NAGLCharges, using NAGL model "
+                                        f"{new_key.extras['partial_charge_method']}, "
+                                        f"applied to topology atom index {topology_atom_index}",
+                                    )
 
                                 elif new_key.extras["handler"] == "preset":
                                     logger.info(

pyproject.toml

@@ -59,7 +59,7 @@ ignore_missing_imports = true
 markers = [
     "slow: marks tests as slow (deselect with '-m \"not slow\"')",
 ]
-addopts = "--cov=openff/interchange --cov-report=xml"
+#addopts = "--cov=openff/interchange --cov-report=xml"
 
 [tool.interrogate]
 ignore-init-method = true