From a7765ee0d980ef7d6f1b3b5ac04923f7fdc1e411 Mon Sep 17 00:00:00 2001
From: Rahil Doshi <rahil.doshi@fau.de>
Date: Tue, 25 Mar 2025 16:05:29 +0100
Subject: [PATCH] Add models for specific enthalpy

---
 src/pymatlib/core/models.py      | 44 +++++++++++++++++++++++++++-----
 src/pymatlib/core/yaml_parser.py | 19 ++++++++++++++
 2 files changed, 56 insertions(+), 7 deletions(-)

diff --git a/src/pymatlib/core/models.py b/src/pymatlib/core/models.py
index 5a3701a..afdbb02 100644
--- a/src/pymatlib/core/models.py
+++ b/src/pymatlib/core/models.py
@@ -75,15 +75,17 @@ def _prepare_material_expressions(*properties):
     """Prepare expressions and collect assignments from material properties."""
     sub_assignments = []
     expressions = []
-
     for prop in properties:
         if isinstance(prop, MaterialProperty):
             sub_assignments.extend(prop.assignments)
             expressions.append(prop.expr)
         else:
             expressions.append(sympy_wrapper(prop))
-
-    return expressions, sub_assignments
+    # If there's only one expression, return it directly instead of in a list
+    if len(expressions) == 1:
+        return expressions[0], sub_assignments
+    else:
+        return expressions, sub_assignments
 
 
 def density_by_thermal_expansion(
@@ -119,10 +121,12 @@ def density_by_thermal_expansion(
         # raise ValueError(f"Thermal expansion coefficient must be between -3e-5 and 0.001, got {thermal_expansion_coefficient}")
 
     try:
-        tec_expr = thermal_expansion_coefficient.expr if isinstance(thermal_expansion_coefficient, MaterialProperty) else sympy_wrapper(thermal_expansion_coefficient)
-        sub_assignments = thermal_expansion_coefficient.assignments if isinstance(thermal_expansion_coefficient, MaterialProperty) else []
-        density = density_base * (1 + tec_expr * (temperature - temperature_base)) ** (-3)
-        return MaterialProperty(density, sub_assignments)
+        tec_expr, sub_assignments \
+            = _prepare_material_expressions(thermal_expansion_coefficient)
+        # tec_expr = thermal_expansion_coefficient.expr if isinstance(thermal_expansion_coefficient, MaterialProperty) else sympy_wrapper(thermal_expansion_coefficient)
+        # sub_assignments = thermal_expansion_coefficient.assignments if isinstance(thermal_expansion_coefficient, MaterialProperty) else []
+        density_expr = density_base * (1 + tec_expr * (temperature - temperature_base)) ** (-3)
+        return MaterialProperty(density_expr, sub_assignments)
     except ZeroDivisionError:
         raise ValueError("Division by zero encountered in density calculation")
 
@@ -157,6 +161,32 @@ def thermal_diffusivity_by_heat_conductivity(
         raise ValueError("Division by zero encountered in thermal diffusivity calculation")
 
 
+# TODO: Add models for specific_enthalpy (h)
+def specific_enthalpy_sensible(
+        temperature: Union[float, sp.Symbol],
+        heat_capacity: Union[float, MaterialProperty]) \
+        -> MaterialProperty:
+
+    heat_capacity_expr, sub_assignments \
+        = _prepare_material_expressions(heat_capacity)
+
+    specific_enthalpy_expr = temperature * heat_capacity_expr
+    return MaterialProperty(specific_enthalpy_expr, sub_assignments)
+
+
+def specific_enthalpy_with_latent_heat(
+        temperature: Union[float, sp.Symbol],
+        heat_capacity: Union[float, MaterialProperty],
+        latent_heat: Union[float, MaterialProperty]) \
+        -> MaterialProperty:
+
+    (heat_capacity_expr, latent_heat_expr), sub_assignments \
+        = _prepare_material_expressions(heat_capacity, latent_heat)
+
+    specific_enthalpy_expr = temperature * heat_capacity_expr + latent_heat_expr
+    return MaterialProperty(specific_enthalpy_expr, sub_assignments)
+
+
 def energy_density(
         density: Union[float, MaterialProperty],
         specific_enthalpy: Union[float, MaterialProperty]) \
diff --git a/src/pymatlib/core/yaml_parser.py b/src/pymatlib/core/yaml_parser.py
index f8974d5..f137b3c 100644
--- a/src/pymatlib/core/yaml_parser.py
+++ b/src/pymatlib/core/yaml_parser.py
@@ -14,6 +14,10 @@ from ruamel.yaml import YAML, constructor, scanner
 from difflib import get_close_matches
 from enum import Enum, auto
 
+from pymatlib.core.models import specific_enthalpy_sensible
+
+from pymatlib.core.models import specific_enthalpy_with_latent_heat
+
 
 class PropertyType(Enum):
     CONSTANT = auto()
@@ -1023,6 +1027,17 @@ class MaterialConfigParser:
                     alloy.heat_capacity
                 )
             },
+            'specific_enthalpy': {
+                'default': lambda: specific_enthalpy_sensible(
+                    T,
+                    alloy.heat_capacity
+                ),
+                'latent_heat_based': lambda: specific_enthalpy_with_latent_heat(
+                    T,
+                    alloy.heat_capacity,
+                    alloy.latent_heat_of_fusion
+                )
+            },
             'energy_density': {
                 'default': lambda: energy_density(
                     alloy.density,
@@ -1046,6 +1061,10 @@ class MaterialConfigParser:
             'thermal_diffusivity': {
                 'default': ['heat_conductivity', 'density', 'heat_capacity'],
             },
+            'specific_enthalpy': {
+                'default': ['heat_capacity'],
+                'latent_heat_based': ['heat_capacity', 'latent_heat_of_fusion'],
+            },
             'energy_density': {
                 'default': ['density', 'specific_enthalpy'],
                 # 'enthalpy_based': ['density', 'specific_enthalpy', 'latent_heat_of_fusion'],
-- 
GitLab