Merge pull request #338 from equinor/option-to-import-component-properties-from-chemicals-database2

option to import component properties from chemicals database

Author: EvenSol

src/neqsim/thermo/thermoTools.py

@@ -265,8 +265,10 @@
 
 """
 
+import importlib
 import logging
-from typing import List, Optional, Union
+from dataclasses import dataclass
+from typing import List, Optional, Tuple, Union
 import jpype
 import pandas
 from jpype.types import *
@@ -315,6 +317,270 @@
 }
 
 
+class ExtendedDatabaseError(Exception):
+    """Raised when a component cannot be resolved in the extended database."""
+
+
+@dataclass
+class _ChemicalComponentData:
+    name: str
+    CAS: str
+    tc: float
+    pc: float
+    omega: float
+    molar_mass: Optional[float] = None
+    normal_boiling_point: Optional[float] = None
+    triple_point_temperature: Optional[float] = None
+    critical_volume: Optional[float] = None
+    critical_compressibility: Optional[float] = None
+
+
+def _create_extended_database_provider():
+    """Create a chemicals database provider."""
+
+    return _ChemicalsDatabaseProvider()
+
+
+class _ChemicalsDatabaseProvider:
+    """Lookup component data from the `chemicals` package."""
+
+    def __init__(self):
+        try:
+            from chemicals.identifiers import CAS_from_any
+        except ImportError as exc:  # pragma: no cover - import guard
+            raise ModuleNotFoundError(
+                "The 'chemicals' package is required to use the extended component database."
+            ) from exc
+
+        self._cas_from_any = CAS_from_any
+        critical = importlib.import_module("chemicals.critical")
+        try:
+            phase_change = importlib.import_module("chemicals.phase_change")
+        except ImportError:  # pragma: no cover - optional submodule
+            phase_change = None
+        try:
+            elements = importlib.import_module("chemicals.elements")
+        except ImportError:  # pragma: no cover - optional submodule
+            elements = None
+
+        self._tc = getattr(critical, "Tc")
+        self._pc = getattr(critical, "Pc")
+        self._omega = getattr(critical, "omega")
+        self._vc = getattr(critical, "Vc", None)
+        self._zc = getattr(critical, "Zc", None)
+        triple_point_candidates = [
+            getattr(critical, "Ttriple", None),
+            getattr(critical, "Tt", None),
+        ]
+        if phase_change is not None:
+            triple_point_candidates.append(getattr(phase_change, "Tt", None))
+        self._triple_point = next(
+            (func for func in triple_point_candidates if func), None
+        )
+        self._tb = (
+            getattr(phase_change, "Tb", None) if phase_change is not None else None
+        )
+        self._molecular_weight = (
+            getattr(elements, "molecular_weight", None)
+            if elements is not None
+            else None
+        )
+
+    def get_component(self, name: str) -> _ChemicalComponentData:
+        cas = self._cas_from_any(name)
+        if not cas:
+            raise ExtendedDatabaseError(
+                f"Component '{name}' was not found in the chemicals database."
+            )
+
+        tc = self._tc(cas)
+        pc = self._pc(cas)
+        omega = self._omega(cas)
+
+        if None in (tc, pc, omega):
+            raise ExtendedDatabaseError(
+                f"Incomplete property data for '{name}' (CAS {cas})."
+            )
+
+        molar_mass = self._call_optional(self._molecular_weight, cas)
+        if molar_mass is not None:
+            molar_mass = float(molar_mass) / 1000.0
+
+        normal_boiling_point = self._call_optional(self._tb, cas)
+        triple_point_temperature = self._call_optional(self._triple_point, cas)
+
+        critical_volume = self._call_optional(self._vc, cas)
+        if critical_volume is not None:
+            critical_volume = float(critical_volume) * 1.0e6  # m^3/mol -> cm^3/mol
+
+        critical_compressibility = self._call_optional(self._zc, cas)
+
+        return _ChemicalComponentData(
+            name=name,
+            CAS=cas,
+            tc=float(tc),
+            pc=float(pc) / 1.0e5,  # chemicals returns pressure in Pa
+            omega=float(omega),
+            molar_mass=molar_mass,
+            normal_boiling_point=(
+                float(normal_boiling_point)
+                if normal_boiling_point is not None
+                else None
+            ),
+            triple_point_temperature=(
+                float(triple_point_temperature)
+                if triple_point_temperature is not None
+                else None
+            ),
+            critical_volume=critical_volume,
+            critical_compressibility=(
+                float(critical_compressibility)
+                if critical_compressibility is not None
+                else None
+            ),
+        )
+
+    @staticmethod
+    def _call_optional(func, cas):
+        if func is None:
+            return None
+        for call in (
+            lambda: func(cas),
+            lambda: func(CASRN=cas),
+        ):
+            try:
+                value = call()
+            except TypeError:
+                continue
+            except Exception:  # pragma: no cover - defensive fallback
+                return None
+            else:
+                return value
+        return None
+
+
+def _get_extended_provider(system):
+    provider = getattr(system, "_extended_database_provider", None)
+    if provider is None:
+        provider = _create_extended_database_provider()
+        system._extended_database_provider = provider  # type: ignore[attr-defined]
+    return provider
+
+
+def _apply_extended_properties(
+    system, component_names: Tuple[str, ...], data: _ChemicalComponentData
+):
+    setter_map = {
+        "CAS": "setCASnumber",
+        "molar_mass": "setMolarMass",
+        "normal_boiling_point": "setNormalBoilingPoint",
+        "triple_point_temperature": "setTriplePointTemperature",
+        "critical_volume": "setCriticalVolume",
+        "critical_compressibility": "setCriticalCompressibilityFactor",
+    }
+
+    for phase_index in range(system.getNumberOfPhases()):
+        try:
+            phase = system.getPhase(phase_index)
+        except Exception:  # pragma: no cover - defensive fallback
+            continue
+        if not hasattr(phase, "hasComponent"):
+            continue
+        component = None
+        for name in component_names:
+            if phase.hasComponent(name):
+                component = phase.getComponent(name)
+                break
+        if component is None:
+            continue
+        for field, setter_name in setter_map.items():
+            value = getattr(data, field, None)
+            if value is None:
+                continue
+            setter = getattr(component, setter_name, None)
+            if setter is None:
+                continue
+            setter(value)
+
+
+def _system_interface_class():
+    """Return the JPype proxy for ``neqsim.thermo.system.SystemInterface``."""
+
+    if not hasattr(_system_interface_class, "_cached"):
+        _system_interface_class._cached = jpype.JClass(  # type: ignore[attr-defined]
+            "neqsim.thermo.system.SystemInterface"
+        )
+    return _system_interface_class._cached  # type: ignore[attr-defined]
+
+
+def _resolve_alias(name: str) -> str:
+    try:
+        return jneqsim.thermo.component.Component.getComponentNameFromAlias(name)
+    except Exception:  # pragma: no cover - defensive alias resolution
+        return name
+
+
+def _has_component_in_database(name: str) -> bool:
+    database = jneqsim.util.database.NeqSimDataBase
+    return database.hasComponent(name) or database.hasTempComponent(name)
+
+
+def _args_look_like_component_properties(args: Tuple[object, ...]) -> bool:
+    return len(args) == 3 and all(isinstance(value, (int, float)) for value in args)
+
+
+@jpype.JImplementationFor("neqsim.thermo.system.SystemInterface")
+class _SystemInterface:
+    def useExtendedDatabase(self, enable: bool = True):
+        """Enable or disable usage of the chemicals based component database."""
+
+        if enable:
+            provider = _create_extended_database_provider()
+            self._use_extended_database = True  # type: ignore[attr-defined]
+            self._extended_database_provider = provider  # type: ignore[attr-defined]
+        else:
+            self._use_extended_database = False  # type: ignore[attr-defined]
+            if hasattr(self, "_extended_database_provider"):
+                delattr(self, "_extended_database_provider")
+        return self
+
+    def addComponent(self, name, amount, *args):  # noqa: N802 - Java signature
+        alias_name = _resolve_alias(name)
+        component_data = None
+
+        if getattr(
+            self, "_use_extended_database", False
+        ) and not _has_component_in_database(alias_name):
+            try:
+                provider = _get_extended_provider(self)
+                component_data = provider.get_component(name)
+            except (ExtendedDatabaseError, ModuleNotFoundError):
+                component_data = None
+
+            if component_data is not None and not _args_look_like_component_properties(
+                args
+            ):
+                if args:
+                    raise NotImplementedError(
+                        "Extended database currently supports components specified in moles (unit='no') "
+                        "without explicit phase targeting or alternative units."
+                    )
+                result = _system_interface_class().addComponent(
+                    self,
+                    name,
+                    float(amount),
+                    component_data.tc,
+                    component_data.pc,
+                    component_data.omega,
+                )
+
+                _apply_extended_properties(self, (alias_name, name), component_data)
+
+                return result
+
+        return _system_interface_class().addComponent(self, name, amount, *args)
+
+
 def fluid(name="srk", temperature=298.15, pressure=1.01325):
     """
     Create a thermodynamic fluid system.
@@ -1100,6 +1366,31 @@ def addComponent(thermoSystem, name, moles, unit="no", phase=-10):
     Returns:
     None
     """
+    alias_name = _resolve_alias(name)
+
+    if getattr(
+        thermoSystem, "_use_extended_database", False
+    ) and not _has_component_in_database(alias_name):
+        try:
+            provider = _get_extended_provider(thermoSystem)
+            component_data = provider.get_component(name)
+        except (ExtendedDatabaseError, ModuleNotFoundError):
+            component_data = None
+        if component_data is not None:
+            if unit != "no" or phase != -10:
+                raise NotImplementedError(
+                    "Extended database currently supports components specified in moles (unit='no') "
+                    "without explicit phase targeting."
+                )
+            thermoSystem.addComponent(
+                name,
+                moles,
+                component_data.tc,
+                component_data.pc,
+                component_data.omega,
+            )
+            return
+
     if phase == -10 and unit == "no":
         thermoSystem.addComponent(name, moles)
     elif phase == -10:

tests/test_extended_database.py

@@ -0,0 +1,77 @@
+import pytest
+
+
+chemicals = pytest.importorskip("chemicals")
+
+import chemicals.critical as critical_data  # type: ignore  # noqa: E402
+from chemicals.critical import Pc, Tc, Vc, Zc, omega  # type: ignore  # noqa: E402
+from chemicals.elements import molecular_weight  # type: ignore  # noqa: E402
+from chemicals.phase_change import Tb  # type: ignore  # noqa: E402
+from chemicals.identifiers import CAS_from_any  # type: ignore  # noqa: E402
+
+from neqsim.thermo.thermoTools import addComponent, fluid
+
+
+def test_use_extended_database_allows_missing_component():
+    system = fluid("srk")
+
+    with pytest.raises(Exception):
+        system.addComponent("dimethylsulfoxide", 1.0)
+
+    system.useExtendedDatabase(True)
+    system.addComponent("dimethylsulfoxide", 1.0)
+
+    component = system.getPhase(0).getComponent("dimethylsulfoxide")
+    cas = CAS_from_any("dimethylsulfoxide")
+
+    assert pytest.approx(component.getTC(), rel=1e-6) == Tc(cas)
+    assert pytest.approx(component.getPC(), rel=1e-6) == Pc(cas) / 1.0e5
+    assert pytest.approx(component.getAcentricFactor(), rel=1e-6) == omega(cas)
+
+    molar_mass = molecular_weight(CASRN=cas)
+    assert molar_mass is not None
+    assert pytest.approx(component.getMolarMass(), rel=1e-6) == molar_mass / 1000.0
+
+    normal_boiling_point = Tb(cas)
+    if normal_boiling_point is not None:
+        assert (
+            pytest.approx(component.getNormalBoilingPoint(), rel=1e-6)
+            == normal_boiling_point
+        )
+
+    critical_volume = Vc(cas)
+    if critical_volume is not None:
+        assert (
+            pytest.approx(component.getCriticalVolume(), rel=1e-6)
+            == critical_volume * 1.0e6
+        )
+
+    critical_compressibility = Zc(cas)
+    if critical_compressibility is not None:
+        assert (
+            pytest.approx(component.getCriticalCompressibilityFactor(), rel=1e-6)
+            == critical_compressibility
+        )
+
+    triple_point_func = getattr(critical_data, "Ttriple", None) or getattr(
+        critical_data, "Tt", None
+    )
+    if triple_point_func is not None:
+        triple_point_temperature = triple_point_func(cas)
+        if triple_point_temperature is not None:
+            assert (
+                pytest.approx(component.getTriplePointTemperature(), rel=1e-6)
+                == triple_point_temperature
+            )
+
+
+def test_module_add_component_uses_extended_database():
+    system = fluid("srk")
+
+    with pytest.raises(Exception):
+        addComponent(system, "dimethylsulfoxide", 1.0)
+
+    system.useExtendedDatabase(True)
+    addComponent(system, "dimethylsulfoxide", 1.0)
+
+    assert system.getPhase(0).hasComponent("dimethylsulfoxide")