Storing meta data in pmd structure (#496)

* storing meta data in pmd structure

* Handle case where forcefield is missing non-bonded force section

* add unit test

* change pytest syntax to capture warning

Author: daico007

foyer/forcefield.py

@@ -29,6 +29,7 @@
     RBTorsionGenerator,
     _convertParameterToNumber,
 )
+from parmed.gromacs.gromacstop import _Defaults
 from pkg_resources import iter_entry_points, resource_filename
 
 import foyer.element as custom_elem
@@ -782,7 +783,7 @@ def apply(
 
         self._apply_typemap(structure, typemap)
 
-        return self.parametrize_system(
+        structure = self.parametrize_system(
             structure=structure,
             references_file=references_file,
             assert_bond_params=assert_bond_params,
@@ -794,6 +795,31 @@ def apply(
             **kwargs,
         )
 
+        # Start storing scaling factors and combining rule to parmed structure
+        # Utilizing parmed gromactop's _Default class
+
+        # Note: nb_func = 1 (LJ) or 2 (Buckingham), foyer forcefield only support 1 at the moment
+        # Combining rule follow GROMACS scheme, where "lorentz" is 2 and "geometric" is 3
+        combining_rules = {"lorentz": 2, "geometric": 3}
+        gen_pairs = "yes" if structure.adjusts else "no"
+        try:
+            lj14scale = self.lj14scale
+            coulomb14scale = self.coulomb14scale
+            structure.defaults = _Defaults(
+                nbfunc=1,
+                comb_rule=combining_rules[self.combining_rule],
+                gen_pairs=gen_pairs,
+                fudgeLJ=lj14scale,
+                fudgeQQ=coulomb14scale,
+            )
+        except AttributeError:
+            warnings.warn(
+                "Missing lj14scale or coulomb14scale, could not set structure metadata."
+            )
+            structure.defaults = None
+
+        return structure
+
     def run_atomtyping(self, structure, use_residue_map=True, **kwargs):
         """Atomtype the topology.
 

foyer/tests/test_forcefield.py

@@ -7,6 +7,7 @@
 import parmed as pmd
 import pytest
 from lxml import etree as ET
+from parmed.gromacs.gromacstop import _Defaults
 from pkg_resources import resource_filename
 
 from foyer import Forcefield, forcefields
@@ -94,6 +95,17 @@ def test_from_parmed(self, oplsaa):
 
         assert ethane.box_vectors == mol2.box_vectors
 
+    def test_structure_meta(self, oplsaa):
+        mol2 = pmd.load_file(get_fn("ethane.mol2"), structure=True)
+        ethane = oplsaa.apply(mol2)
+
+        assert isinstance(ethane.defaults, _Defaults)
+        assert ethane.defaults.nbfunc == 1
+        assert ethane.defaults.comb_rule == 3
+        assert ethane.defaults.fudgeLJ == oplsaa.lj14scale
+        assert ethane.defaults.fudgeQQ == oplsaa.coulomb14scale
+        assert ethane.defaults.gen_pairs == "yes"
+
     @pytest.mark.skipif(not has_mbuild, reason="mbuild is not installed")
     def test_from_mbuild(self, oplsaa):
         import mbuild as mb
@@ -426,7 +438,10 @@ def test_allow_empty_def(self):
         ethane = mb.load(get_fn("ethane.mol2"))
         with pytest.warns(ValidationWarning):
             ff = Forcefield(forcefield_files=get_fn("empty_def.xml"))
-        ff.apply(ethane)
+
+        with pytest.warns(UserWarning):
+            typed = ff.apply(ethane)
+        assert typed.defaults is None
 
     @pytest.mark.skipif(not has_mbuild, reason="mbuild is not installed")
     def test_assert_bonds(self):