[PH-Flash] Add mixture enthalpy models with cited caloric data and baseline flash tests

CMakeLists_files.cmake

@@ -709,7 +709,9 @@ list(APPEND DUNE_TEST_SOURCE_FILES
   tests/material/test_ncpflash.cpp
   tests/material/test_pengrobinson.cpp
   tests/material/test_tabulation.cpp
+  tests/material/test_mixture_enthalpy.cpp
   tests/material/test_threecomponents_ptflash.cpp
+  tests/material/test_twocomponents_ptflash.cpp
 )
 
 if(dune-common_FOUND)
@@ -1312,6 +1314,7 @@ list(APPEND PUBLIC_HEADER_FILES
   opm/material/components/BrineDynamic.hpp
   opm/material/components/C1.hpp
   opm/material/components/C10.hpp
+  opm/material/components/ComponentCp.hpp
   opm/material/components/CO2.hpp
   opm/material/components/CO2Tables.hpp
   opm/material/components/Component.hpp
@@ -1334,8 +1337,10 @@ list(APPEND PUBLIC_HEADER_FILES
   opm/material/components/iapws/Region4.hpp
   opm/material/constraintsolvers/CompositionFromFugacities.hpp
   opm/material/constraintsolvers/ComputeFromReferencePhase.hpp
+  opm/material/constraintsolvers/IdealGasCaloricData.hpp
   opm/material/constraintsolvers/ImmiscibleFlash.hpp
   opm/material/constraintsolvers/MiscibleMultiPhaseComposition.hpp
+  opm/material/constraintsolvers/MixtureEnthalpy.hpp
   opm/material/constraintsolvers/NcpFlash.hpp
   opm/material/constraintsolvers/PTFlash.hpp
   opm/material/densead/DynamicEvaluation.hpp

opm/material/components/C1.hpp

@@ -28,6 +28,7 @@
 #define OPM_C1_HPP
 
 #include "Component.hpp"
+#include "ComponentCp.hpp"
 
 #include <opm/material/IdealGas.hpp>
 #include <opm/material/common/MathToolbox.hpp>
@@ -42,7 +43,11 @@ namespace Opm
 /*!
  * \ingroup Components
  *
- * \brief Properties of pure molecular methane \f$C_1\f$.
+ * \brief Properties of pure molecular methane \f$C_1\f$
+ *        (CAS 74-82-8, formula CH\f$_4\f$).
+ *
+ * Constants cross-checked against the CoolProp reference fluid "Methane"
+ * (reference EoS: Setzmann & Wagner, J. Phys. Chem. Ref. Data 20 (1991) 1061).
  *
  * \tparam Scalar The type used for scalar values
  */
@@ -53,7 +58,7 @@ class C1 : public Component<Scalar, C1<Scalar> >
 
 public:
     /*!
-     * \brief A human readable name for NDecane.
+     * \brief A human readable name for methane.
      */
     static std::string_view name()
     { return "C1"; }
@@ -86,9 +91,34 @@ class C1 : public Component<Scalar, C1<Scalar> >
      */
     static Scalar acentricFactor() { return 0.011; }
 
+    /*!
+     * \brief Returns the temperature \f$\mathrm{[K]}\f$ at the triple point of methane.
+     */
+    static Scalar tripleTemperature()
+    { return 90.6941; /* [K] — Setzmann & Wagner (1991), via CoolProp 8.0.0 */ }
 
+    /*!
+     * \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at the triple point of methane.
+     */
+    static Scalar triplePressure()
+    { return 1.1696e4; /* [Pa] — Setzmann & Wagner (1991), via CoolProp 8.0.0 */ }
 
-
+    /*!
+     * \brief Cubic ideal-gas heat-capacity polynomial of methane \f$\mathrm{[J/(mol\ K)]}\f$.
+     *
+     * Least-squares fit to the ideal-gas part of the reference EoS
+     * (Setzmann & Wagner 1991), window 250–600 K, RMS 0.024 /
+     * max 0.063 J/(mol K); sampled via CoolProp 8.0.0 (extraction tool
+     * only). Outside the window the cubic extrapolates — refit rather than
+     * trust it there.
+     *
+     * This is the CALORIC (ideal-gas) identity consumed by the compositional
+     * mixture-enthalpy model; it deliberately does NOT implement the
+     * Component<> gasEnthalpy/gasHeatCapacity slots (those are real-fluid
+     * correlations where implemented).
+     */
+    static constexpr ComponentCp<Scalar> idealGasHeatCapacityPolynomial()
+    { return {40.1503, -8.47372e-2, 2.93012e-4, -1.96125e-7}; }
 };
 
 } // namespace Opm

opm/material/components/C10.hpp

@@ -28,6 +28,7 @@
 #define OPM_C10_HPP
 
 #include "Component.hpp"
+#include "ComponentCp.hpp"
 
 #include <opm/material/IdealGas.hpp>
 #include <opm/material/common/MathToolbox.hpp>
@@ -41,7 +42,11 @@ namespace Opm
 /*!
  * \ingroup Components
  *
- * \brief Properties of pure molecular n-Decane \f$C_10\f$.
+ * \brief Properties of pure molecular n-Decane \f$C_10\f$
+ *        (CAS 124-18-5, formula C\f$_{10}\f$H\f$_{22}\f$).
+ *
+ * Constants cross-checked against the CoolProp reference fluid "n-Decane"
+ * (reference EoS: Lemmon & Span, J. Chem. Eng. Data 51 (2006) 785).
  *
  * \tparam Scalar The type used for scalar values
  */
@@ -85,7 +90,34 @@ class C10 : public Component<Scalar, C10<Scalar> >
      */
     static Scalar acentricFactor() { return 0.488; }
 
+    /*!
+     * \brief Returns the temperature \f$\mathrm{[K]}\f$ at the triple point of n-Decane.
+     */
+    static Scalar tripleTemperature()
+    { return 243.5; /* [K] — Lemmon & Span (2006), via CoolProp 8.0.0 */ }
+
+    /*!
+     * \brief Returns the pressure \f$\mathrm{[Pa]}\f$ at the triple point of n-Decane.
+     */
+    static Scalar triplePressure()
+    { return 1.404; /* [Pa] — Lemmon & Span (2006), via CoolProp 8.0.0 */ }
 
+    /*!
+     * \brief Cubic ideal-gas heat-capacity polynomial of n-Decane \f$\mathrm{[J/(mol\ K)]}\f$.
+     *
+     * Least-squares fit to the ideal-gas part of the reference EoS
+     * (Lemmon & Span 2006), window 250–600 K, RMS 0.266 /
+     * max 0.792 J/(mol K); sampled via CoolProp 8.0.0 (extraction tool
+     * only). Outside the window the cubic extrapolates — refit rather than
+     * trust it there.
+     *
+     * This is the CALORIC (ideal-gas) identity consumed by the compositional
+     * mixture-enthalpy model; it deliberately does NOT implement the
+     * Component<> gasEnthalpy/gasHeatCapacity slots (those are real-fluid
+     * correlations where implemented).
+     */
+    static constexpr ComponentCp<Scalar> idealGasHeatCapacityPolynomial()
+    { return {79.4791, 3.1066e-1, 9.88317e-4, -1.00245e-6}; }
 };
 
 } // namespace Opm

opm/material/components/ComponentCp.hpp

@@ -0,0 +1,82 @@
+// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
+// vi: set et ts=4 sw=4 sts=4:
+/*
+  Copyright 2026 Equinor ASA.
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM is free software: you can redistribute it and/or modify
+  it under the terms of the GNU General Public License as published by
+  the Free Software Foundation, either version 2 of the License, or
+  (at your option) any later version.
+
+  OPM is distributed in the hope that it will be useful,
+  but WITHOUT ANY WARRANTY; without even the implied warranty of
+  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+  GNU General Public License for more details.
+
+  You should have received a copy of the GNU General Public License
+  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
+
+  Consult the COPYING file in the top-level source directory of this
+  module for the precise wording of the license and the list of
+  copyright holders.
+*/
+/*!
+ * \file
+ *
+ * \brief The cubic ideal-gas heat-capacity polynomial and its closed-form
+ *        enthalpy integral — the caloric EQUATIONS, species-blind.
+ *
+ * The coefficients (a species' caloric identity) live on the component
+ * classes (e.g. C1::idealGasHeatCapacityPolynomial()); the name-keyed lookup
+ * and the enthalpy reference datum live in IdealGasCaloricData. This header
+ * holds only the mathematics shared by all of them.
+ *
+ * Units are SI throughout: temperature [K], molar heat capacity [J/(mol K)],
+ * molar enthalpy [J/mol].
+ */
+#ifndef OPM_COMPONENT_CP_HPP
+#define OPM_COMPONENT_CP_HPP
+
+#include <array>
+
+namespace Opm {
+
+/*!
+ * \brief Cubic ideal-gas heat-capacity polynomial of one component:
+ *        cp(T) = c0 + c1*T + c2*T^2 + c3*T^3   [J/(mol K)]
+ */
+template <class Scalar>
+struct ComponentCp {
+    Scalar c0, c1, c2, c3;
+
+    //! cp(T) [J/(mol K)]. Generic in the evaluation type (double or AD).
+    template <class Eval>
+    Eval heatCapacity(const Eval& T) const
+    {
+        return c0 + c1*T + c2*T*T + c3*T*T*T;
+    }
+
+    /*!
+     * \brief Ideal-gas enthalpy h(T) = int_{T0}^{T} cp dT' [J/mol],
+     *        in closed form. h(T0) = 0 by construction.
+     */
+    template <class Eval>
+    Eval enthalpyIntegral(const Eval& T, const Scalar T0) const
+    {
+        return c0*(T - T0)
+             + c1/2*(T*T - T0*T0)
+             + c2/3*(T*T*T - T0*T0*T0)
+             + c3/4*(T*T*T*T - T0*T0*T0*T0);
+    }
+};
+
+//! Per-component cp table for an N-component fluid system, indexed like the
+//! fluid system's component indices.
+template <class Scalar, int numComponents>
+using CpTable = std::array<ComponentCp<Scalar>, numComponents>;
+
+} // namespace Opm
+
+#endif // OPM_COMPONENT_CP_HPP

opm/material/components/SimpleCO2.hpp

@@ -32,6 +32,7 @@
 
 #include <opm/material/IdealGas.hpp>
 #include <opm/material/components/Component.hpp>
+#include <opm/material/components/ComponentCp.hpp>
 
 #include <cmath>
 #include <string_view>
@@ -42,6 +43,11 @@ namespace Opm {
  * \ingroup Components
  *
  * \brief A simplistic class representing the \f$CO_2\f$ fluid properties
+ *        (CAS 124-38-9, formula CO\f$_2\f$).
+ *
+ * Constants cross-checked against the CoolProp reference fluid
+ * "CarbonDioxide" (reference EoS: Span & Wagner, J. Phys. Chem. Ref. Data
+ * 25 (1996) 1509).
  *
  * \tparam Scalar The type used for scalar values
  */
@@ -98,6 +104,23 @@ class SimpleCO2 : public Component<Scalar, SimpleCO2<Scalar> >
     static Scalar triplePressure()
     { return 5.11e5; /* [N/m^2] */ }
 
+    /*!
+     * \brief Cubic ideal-gas heat-capacity polynomial of \f$CO_2\f$ \f$\mathrm{[J/(mol\ K)]}\f$.
+     *
+     * Least-squares fit to the ideal-gas part of the reference EoS
+     * (Span & Wagner 1996), window 250–600 K, RMS 0.003 /
+     * max 0.014 J/(mol K); sampled via CoolProp 8.0.0 (extraction tool
+     * only). Outside the window the cubic extrapolates — refit rather than
+     * trust it there.
+     *
+     * This is the CALORIC (ideal-gas) identity consumed by the compositional
+     * mixture-enthalpy model; it deliberately does NOT implement the
+     * Component<> gasEnthalpy/gasHeatCapacity slots (those are real-fluid
+     * correlations where implemented).
+     */
+    static constexpr ComponentCp<Scalar> idealGasHeatCapacityPolynomial()
+    { return {18.2687, 8.36359e-2, -7.75148e-5, 3.14088e-8}; }
+
     /*!
      * \copydoc Component::gasIsCompressible
      */