added interaction groupings to gromacs parser

atomisticparsers/gromacs/parser.py

@@ -833,10 +833,53 @@ def parse_method(self):
             self.logger.error('Error parsing interactions.')
 
         interactions = self.traj_parser.get_interactions()
-        for interaction in interactions:
+        #  TODO The section below is duplicated in the lammps and gromacs parsers. We should really move them
+        #       to the MDAnalysis parser, but you must be careful because there are other contributions to interactions
+        #       from other (sub)parsers
+        interaction_key_list = Interaction.__dict__.keys()
+        interaction_dict = {}
+        interaction_keys_remove = ['__module__', '__doc__', 'm_def']
+        interaction_key_list = [key for key in interaction_key_list if key not in interaction_keys_remove]
+        for interaction_key in interaction_key_list:
+            interaction_dict[interaction_key] = np.array([interaction.get(interaction_key) for interaction in interactions], dtype=object)
+        interaction_dict = {key: val for key, val in interaction_dict.items()}
+        interaction_types = np.unique(interaction_dict['type']) if interaction_dict.get('type') is not None else []
+        for interaction_type in interaction_types:
             sec_interaction = sec_model.m_create(Interaction)
-            for key, val in interaction.items():
-                setattr(sec_interaction, key, val)
+            interaction_indices = np.where(interaction_dict['type'] == interaction_type)[0]
+            sec_interaction.type = interaction_type
+            sec_interaction.n_inter = len(interaction_indices)
+            sec_interaction.n_atoms
+            n_atoms = interaction_dict.get('n_atoms')[interaction_indices][0]  # assuming the number of atoms per interaction type is fixed!
+            for key, val in interaction_dict.items():
+                if key == 'type':
+                    continue
+                interaction_vals = val[interaction_indices]
+                if type(interaction_vals[0]).__name__ == 'ndarray':
+                    interaction_vals = np.array([vals.tolist() for vals in interaction_vals], dtype=object)
+                if interaction_vals.all() is None:
+                    continue
+                if key == 'parameters':
+                    interaction_vals = interaction_vals.tolist()
+                elif key == 'n_atoms':
+                    interaction_vals = interaction_vals[0]
+                try:
+                    setattr(sec_interaction, key, interaction_vals)
+                except Exception:
+                    self.logger.warning(
+                        'Some issue trying to store ' + key + 'in Interactions section.'
+                        ' Possibly a data type problem. Ignoring these values.')
+
+            if not sec_interaction.get('n_atoms'):
+                sec_interaction.n_atoms = len(sec_interaction.get('atom_indices')[0]) if sec_interaction.get('atom_indices') is not None else None
+
+        # OLD VERSION WITHOUT GROUPINGS
+        # for interaction in interactions:
+        #     sec_interaction = sec_model.m_create(Interaction)
+        #     for key, val in interaction.items():
+        #         print(key)
+        #         print(val)
+        #         setattr(sec_interaction, key, val)
 
         input_parameters = self.log_parser.get('input_parameters', {})
         sec_force_calculations = sec_force_field.m_create(ForceCalculations)

tests/test_gromacsparser.py

@@ -103,9 +103,16 @@ def test_md_verbose(parser):
 
     sec_method = sec_run.method
     assert len(sec_method) == 1
-    assert len(sec_method[0].force_field.model[0].contributions) == 1127
-    assert sec_method[0].force_field.model[0].contributions[0].type == 'angle'
-    assert sec_method[0].force_field.model[0].contributions[1120].parameters[1] == 575.0
+    # assert len(sec_method[0].force_field.model[0].contributions) == 1127
+    # assert sec_method[0].force_field.model[0].contributions[0].type == 'angle'
+    # assert sec_method[0].force_field.model[0].contributions[1120].parameters[1] == 575.0
+    assert len(sec_method[0].force_field.model[0].contributions) == 8
+    assert sec_method[0].force_field.model[0].contributions[6].type == 'bond'
+    assert sec_method[0].force_field.model[0].contributions[6].n_inter == 1017
+    assert sec_method[0].force_field.model[0].contributions[6].n_atoms == 2
+    assert sec_method[0].force_field.model[0].contributions[6].atom_labels[10][0] == 'C'
+    assert sec_method[0].force_field.model[0].contributions[6].atom_indices[100][1] == 141
+    assert sec_method[0].force_field.model[0].contributions[6].parameters[858] == approx(0.9999996193044006)
     assert sec_method[0].force_field.force_calculations.vdw_cutoff.magnitude == 1.2e-09
     assert sec_method[0].force_field.force_calculations.vdw_cutoff.units == 'meter'
     assert sec_method[0].force_field.force_calculations.coulomb_type == 'particle_mesh_ewald'