Added parsing of , conversion to

Author: Bernadette-Mohr

src/nomad_simulation_parsers/parsers/lammps/parser.py

@@ -125,6 +125,42 @@ def _extract_pressure_params(
             _pressure_params = None
         return _pressure_params
 
+    def _extract_modulus(self, fix_cmd: list[str]) -> np.ndarray | None:
+        """Extract modulus (bulk modulus) from fix command and convert to compressibility.
+
+        LAMMPS uses 'modulus' keyword in fix npt/nph commands to specify bulk modulus.
+        Compressibility = 1 / modulus (in consistent units).
+
+        Args:
+            fix_cmd: Fix command arguments list
+
+        Returns:
+            3x3 compressibility matrix (1/pascal) or None if modulus not found
+        """
+        result = None
+        if 'modulus' not in fix_cmd:
+            return result
+
+        try:
+            modulus_idx = fix_cmd.index('modulus')
+            # Modulus value follows the 'modulus' keyword
+            modulus_value = float(fix_cmd[modulus_idx + 1])
+
+            # Apply pressure units to modulus (same units as pressure in LAMMPS)
+            modulus_with_units = self.apply_unit(modulus_value, 'pressure')
+
+            # Compressibility = 1 / modulus
+            # Units: 1/pascal = 1 / (kg/(m·s²)) = m·s²/kg = Pa⁻¹
+            compressibility_value = 1.0 / modulus_with_units.magnitude
+            compressibility_unit = 1.0 / modulus_with_units.units
+
+            # Create 3x3 diagonal matrix (isotropic compressibility)
+            result = np.eye(3) * compressibility_value * compressibility_unit
+        except (ValueError, IndexError, TypeError, ZeroDivisionError):
+            result = None
+
+        return result
+
     def _extract_timestep(self) -> Any | None:
         """Extract integration timestep from log parser."""
         timestep_value = self._log_parser.get('timestep')
@@ -183,6 +219,51 @@ def _set_parser_masses(self, masses: np.ndarray | None) -> None:
             if isinstance(parser, TrajParser):
                 parser.masses = masses
 
+    def _extract_integrator_type(self) -> str | None:
+        """Extract integrator type from run_style command.
+
+        LAMMPS run_style options: verlet, verlet/split, respa, respa/omp
+        Maps to NOMAD integrator_type enum.
+
+        Returns:
+            Integrator type string or None if not found
+        """
+        run_style = self._log_parser.get('run_style')
+        if run_style is None:
+            return None
+
+        # Handle nested list structure from multiple run_style commands
+        # E.g., [['verlet'], ['respa', '4', '2', ...]] -> take last command
+        if isinstance(run_style, list) and len(run_style) > 0:
+            if isinstance(run_style[0], list):
+                # Nested list: take last command
+                run_style = run_style[-1]
+
+        # Extract style name from command arguments
+        # E.g., ['respa', '4', '2', ...] -> 'respa'
+        #       'verlet' -> 'verlet'
+        if isinstance(run_style, list) and len(run_style) > 0:
+            style = run_style[0]
+        else:
+            style = run_style
+
+        if style is None:
+            return None
+
+        result = None
+        style_lower = str(style).lower()
+
+        # Map LAMMPS run_style to NOMAD integrator_type
+        integrator_map = {
+            'verlet': 'velocity_verlet',
+            'verlet/split': 'velocity_verlet_split',
+            'respa': 'respa',  # reversible reference system propagator algorithm
+            'respa/omp': 'respa',
+        }
+
+        result = integrator_map.get(style_lower, 'velocity_verlet')  # default fallback
+        return result
+
     def _extract_barostat_settings(
         self, fix_commands: list[list[str]] | None
     ) -> BarostatParameters | None:
@@ -228,6 +309,16 @@ def _extract_barostat_settings(
                         pressure_params
                     )
                 break
+            # Extract modulus and convert to compressibility
+            compressibility = self._extract_modulus(fix_cmd)
+            if compressibility is not None:
+                params.compressibility = compressibility
+
+            # Extract modulus and convert to compressibility
+            compressibility = self._extract_modulus(fix_cmd)
+            if compressibility is not None:
+                params.compressibility = compressibility
+
             return params
         return None
 
@@ -328,7 +419,10 @@ def parse_method(self, simulation: Simulation) -> None:
             return
 
         method = MolecularDynamicsMethod()
-        method.integrator_type = 'velocity_verlet'
+
+        integrator_type = self._extract_integrator_type()
+        if integrator_type is not None:
+            method.integrator_type = integrator_type
 
         timestep = self._extract_timestep()
         if timestep is not None:

tests/parsers/test_lammps_parser.py

@@ -26,6 +26,7 @@
 import pytest
 from nomad.client import normalize_all
 from nomad.datamodel import EntryArchive, EntryMetadata
+from nomad.units import ureg
 from nomad.utils import get_logger
 
 from nomad_simulation_parsers.parsers.lammps.file_parsers import LogParser
@@ -836,6 +837,74 @@ def test_md_method_barostat_extraction():
     assert barostat is None
 
 
+def test_md_method_barostat_modulus_extraction():
+    """Test modulus extraction and conversion to compressibility."""
+
+    writer = LammpsArchiveWriter()
+    writer._log_parser = LogParser()
+    # Initialize units system for apply_unit to work
+    writer._log_parser._units = {
+        'pressure': ureg.atmosphere,
+        'time': ureg.femtosecond,
+    }
+
+    # Test NPT with modulus keyword (bulk modulus = 20000 atm)
+    # compressibility should be 1/20000 atm^-1
+    fix_npt_modulus = [
+        [
+            '1',
+            'all',
+            'npt',
+            'temp',
+            '300',
+            '300',
+            '100',
+            'iso',
+            '1.0',
+            '1.0',
+            '1000.0',
+            'modulus',
+            '20000.0',
+        ]
+    ]
+    barostat = writer._extract_barostat_settings(fix_npt_modulus)
+    assert barostat is not None
+    assert barostat.compressibility is not None
+
+    # Check compressibility value (1/20000 atm in SI units)
+    # 1 atm = 101325 Pa, so 1/20000 atm^-1 = 1/(20000*101325) Pa^-1
+    expected_compressibility = 1.0 / (20000.0 * 101325.0)  # Pa^-1
+    assert barostat.compressibility[0, 0].magnitude == pytest.approx(
+        expected_compressibility, rel=1e-6
+    )
+    assert barostat.compressibility[1, 1].magnitude == pytest.approx(
+        expected_compressibility, rel=1e-6
+    )
+    assert barostat.compressibility[2, 2].magnitude == pytest.approx(
+        expected_compressibility, rel=1e-6
+    )
+
+    # Test NPH with modulus
+    fix_nph_modulus = [
+        ['1', 'all', 'nph', 'iso', '1.0', '1.0', '1000.0', 'modulus', '10000.0']
+    ]
+    barostat = writer._extract_barostat_settings(fix_nph_modulus)
+    assert barostat is not None
+    assert barostat.compressibility is not None
+    expected_compressibility = 1.0 / (10000.0 * 101325.0)
+    assert barostat.compressibility[0, 0].magnitude == pytest.approx(
+        expected_compressibility, rel=1e-6
+    )
+
+    # Test without modulus (should be None)
+    fix_npt_no_modulus = [
+        ['1', 'all', 'npt', 'temp', '300', '300', '100', 'iso', '1.0', '1.0', '1000.0']
+    ]
+    barostat = writer._extract_barostat_settings(fix_npt_no_modulus)
+    assert barostat is not None
+    assert barostat.compressibility is None
+
+
 def test_md_method_ensemble_detection():
     """Test ensemble type detection."""
 
@@ -1014,3 +1083,61 @@ def test_md_method_helper_methods():
     writer._log_parser._results = {'thermo': [50, 100, 200]}
     thermo_freq = writer._extract_thermo_frequency()
     assert thermo_freq == 200
+
+
+def test_md_method_integrator_type_extraction():
+    """Test integrator type extraction from run_style command."""
+
+    class MockLogParser:
+        """Mock LogParser for testing."""
+
+        def __init__(self):
+            self._results = {}
+
+        def get(self, key, default=None):
+            return self._results.get(key, default)
+
+    writer = LammpsArchiveWriter()
+    writer._log_parser = MockLogParser()
+
+    # Test verlet
+    writer._log_parser._results = {'run_style': 'verlet'}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'velocity_verlet'
+
+    # Test verlet/split
+    writer._log_parser._results = {'run_style': 'verlet/split'}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'velocity_verlet_split'
+
+    # Test respa
+    writer._log_parser._results = {'run_style': 'respa'}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'respa'
+
+    # Test respa/omp
+    writer._log_parser._results = {'run_style': 'respa/omp'}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'respa'
+
+    # Test with command args (respa with levels)
+    writer._log_parser._results = {
+        'run_style': ['respa', '4', '2', '2', '2', 'bond', '1', 'kspace', '4']
+    }
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'respa'
+
+    # Test with nested list of commands (take last)
+    writer._log_parser._results = {'run_style': [['verlet'], ['respa', '4', '2']]}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'respa'
+
+    # Test None case
+    writer._log_parser._results = {}
+    integrator = writer._extract_integrator_type()
+    assert integrator is None
+
+    # Test unknown style (fallback to velocity_verlet)
+    writer._log_parser._results = {'run_style': 'unknown_style'}
+    integrator = writer._extract_integrator_type()
+    assert integrator == 'velocity_verlet'