[DualPorosity] (1/5) Add the RUNSPEC porosity model

CMakeLists_files.cmake

@@ -90,6 +90,7 @@ list(APPEND MAIN_SOURCE_FILES
   opm/input/eclipse/EclipseState/EclipseState.cpp
   opm/input/eclipse/EclipseState/EndpointScaling.cpp
   opm/input/eclipse/EclipseState/Phase.cpp
+  opm/input/eclipse/EclipseState/PorosityModel.cpp
   opm/input/eclipse/EclipseState/Runspec.cpp
   opm/input/eclipse/EclipseState/TracerConfig.cpp
   opm/input/eclipse/EclipseState/WagHysteresisConfig.cpp
@@ -1005,6 +1006,7 @@ list(APPEND PUBLIC_HEADER_FILES
   opm/input/eclipse/EclipseState/InitConfig/FoamConfig.hpp
   opm/input/eclipse/EclipseState/InitConfig/InitConfig.hpp
   opm/input/eclipse/EclipseState/Phase.hpp
+  opm/input/eclipse/EclipseState/PorosityModel.hpp
   opm/input/eclipse/EclipseState/Runspec.hpp
   opm/input/eclipse/EclipseState/SimulationConfig/BCConfig.hpp
   opm/input/eclipse/EclipseState/SimulationConfig/DatumDepth.hpp

opm/input/eclipse/EclipseState/Grid/FieldProps.hpp

@@ -134,6 +134,7 @@ namespace ALIAS {
 
 namespace GRID {
 static const std::unordered_map<std::string, keyword_info<double>> double_keywords = {{"DISPERC",keyword_info<double>{}.unit_string("Length")},
+                                                                                      {"SIGMAV",  keyword_info<double>{}.unit_string("1/Length*Length")},
                                                                                       {"MINPVV",  keyword_info<double>{}.init(0.0).unit_string("ReservoirVolume").global_kw(true)},
                                                                                       {"MULTPV",  keyword_info<double>{}.init(1.0).mult(true)},
                                                                                       {"NTG",     keyword_info<double>{}.init(1.0)},

opm/input/eclipse/EclipseState/PorosityModel.cpp

@@ -0,0 +1,77 @@
+/*
+  Copyright 2026 TNO
+
+  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 3 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/>.
+*/
+
+#include <opm/input/eclipse/EclipseState/PorosityModel.hpp>
+
+#include <opm/input/eclipse/Deck/Deck.hpp>
+
+#include <opm/input/eclipse/Parser/ParserKeywords/D.hpp>
+#include <opm/input/eclipse/Parser/ParserKeywords/N.hpp>
+
+namespace Opm {
+
+PorosityModel::PorosityModel(const Deck& deck)
+{
+    if (deck.hasKeyword<ParserKeywords::DUALPERM>()) {
+        this->m_type = Type::DualPermeability;
+    } else if (deck.hasKeyword<ParserKeywords::DUALPORO>()) {
+        this->m_type = Type::DualPorosity;
+    }
+
+    if (deck.hasKeyword<ParserKeywords::NODPPM>()) {
+        this->m_scale_fracture_perm = false;
+    }
+}
+
+PorosityModel::Type PorosityModel::type() const noexcept
+{
+    return this->m_type;
+}
+
+bool PorosityModel::dualContinuum() const noexcept
+{
+    return this->m_type != Type::SinglePorosity;
+}
+
+bool PorosityModel::dualPermeability() const noexcept
+{
+    return this->m_type == Type::DualPermeability;
+}
+
+bool PorosityModel::fracturePermeabilityScalingActive() const noexcept
+{
+    return this->dualContinuum() && this->m_scale_fracture_perm;
+}
+
+PorosityModel PorosityModel::serializationTestObject()
+{
+    PorosityModel result;
+    result.m_type = Type::DualPermeability;
+    result.m_scale_fracture_perm = false;
+
+    return result;
+}
+
+bool PorosityModel::operator==(const PorosityModel& other) const
+{
+    return (this->m_type == other.m_type)
+        && (this->m_scale_fracture_perm == other.m_scale_fracture_perm);
+}
+
+} // namespace Opm

opm/input/eclipse/EclipseState/PorosityModel.hpp

@@ -0,0 +1,102 @@
+/*
+  Copyright 2026 TNO
+
+  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 3 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/>.
+*/
+
+#ifndef OPM_POROSITY_MODEL_HPP
+#define OPM_POROSITY_MODEL_HPP
+
+namespace Opm {
+
+class Deck;
+
+/// Porosity model of a simulation run, as selected in the RUNSPEC section.
+///
+/// A run resolves either one continuum per cell (single porosity, the
+/// default) or two: the rock matrix and the fracture system.  DUALPORO
+/// selects the dual-porosity model, in which fluid moves between fracture
+/// cells and between each matrix cell and its fracture twin, and DUALPERM
+/// the dual-permeability model, which adds flow between matrix cells.  In
+/// the dual-continuum models the fracture permeabilities are scaled by the
+/// fracture porosity unless NODPPM switches the scaling off.
+///
+/// This is the selection the EGRID file header records to distinguish
+/// single-porosity, dual-porosity and dual-permeability grids.
+class PorosityModel
+{
+public:
+    /// Which continua the run resolves.
+    enum class Type {
+        /// One continuum per cell.
+        SinglePorosity,
+
+        /// Matrix and fracture continua, with flow between fracture cells
+        /// and between each matrix cell and its fracture twin.
+        DualPorosity,
+
+        /// Dual porosity with flow between matrix cells as well.
+        DualPermeability,
+    };
+
+    /// Single porosity.
+    PorosityModel() = default;
+
+    /// Porosity model selected by the RUNSPEC keywords of a deck.
+    ///
+    /// \param[in] deck Input deck.  DUALPERM takes precedence over DUALPORO.
+    explicit PorosityModel(const Deck& deck);
+
+    /// Selected model.
+    Type type() const noexcept;
+
+    /// Whether the run resolves both a matrix and a fracture continuum,
+    /// i.e., whether the model is dual porosity or dual permeability.
+    bool dualContinuum() const noexcept;
+
+    /// Whether the run resolves flow between matrix cells, i.e., whether
+    /// the model is dual permeability.
+    bool dualPermeability() const noexcept;
+
+    /// Whether the fracture permeabilities are scaled by the fracture
+    /// porosity.  Active in the dual-continuum models unless NODPPM is
+    /// specified, never in single porosity.  Every consumer of the scaling
+    /// rule, well connection factors and cell transmissibilities alike,
+    /// must take the answer from here rather than derive its own.
+    bool fracturePermeabilityScalingActive() const noexcept;
+
+    static PorosityModel serializationTestObject();
+
+    bool operator==(const PorosityModel& other) const;
+
+    template <class Serializer>
+    void serializeOp(Serializer& serializer)
+    {
+        serializer(this->m_type);
+        serializer(this->m_scale_fracture_perm);
+    }
+
+private:
+    Type m_type{Type::SinglePorosity};
+
+    /// Fracture permeability scaling as requested by the deck.  Only
+    /// meaningful in the dual-continuum models; NODPPM turns it off.
+    bool m_scale_fracture_perm{true};
+};
+
+} // namespace Opm
+
+#endif // OPM_POROSITY_MODEL_HPP

opm/input/eclipse/EclipseState/Runspec.cpp

@@ -823,6 +823,7 @@ Runspec::Runspec(const Deck& deck)
     , m_mechsolver (deck)
     , m_tracers    (deck)
     , m_geochem    (deck)
+    , m_porosity_model(deck)
     , m_co2storage (false)
     , m_co2sol     (false)
     , m_h2sol      (false)
@@ -1010,6 +1011,7 @@ Runspec Runspec::serializationTestObject()
     result.m_temp = true;
     result.m_biof = true;
     result.m_geochem = Geochem::serializationTestObject();
+    result.m_porosity_model = PorosityModel::serializationTestObject();
 
     return result;
 }
@@ -1135,6 +1137,11 @@ bool Runspec::biof() const noexcept
     return this->m_biof;
 }
 
+const PorosityModel& Runspec::porosityModel() const noexcept
+{
+    return this->m_porosity_model;
+}
+
 std::time_t Runspec::start_time() const noexcept
 {
     return this->m_start_time;
@@ -1186,6 +1193,7 @@ bool Runspec::rst_cmp(const Runspec& full_spec, const Runspec& rst_spec)
         full_spec.m_temp == rst_spec.m_temp &&
         full_spec.m_biof == rst_spec.m_biof &&
         full_spec.m_geochem == rst_spec.m_geochem &&
+        full_spec.m_porosity_model == rst_spec.m_porosity_model &&
         Welldims::rst_cmp(full_spec.wellDimensions(), rst_spec.wellDimensions());
 }
 
@@ -1218,6 +1226,7 @@ bool Runspec::operator==(const Runspec& data) const
         && (this->m_temp == data.m_temp)
         && (this->m_biof == data.m_biof)
         && (this->m_geochem == data.m_geochem)
+        && (this->m_porosity_model == data.m_porosity_model)
         ;
 }
 

opm/input/eclipse/EclipseState/Runspec.hpp

@@ -23,6 +23,7 @@
 
 #include <opm/input/eclipse/EclipseState/EndpointScaling.hpp>
 #include <opm/input/eclipse/EclipseState/Phase.hpp>
+#include <opm/input/eclipse/EclipseState/PorosityModel.hpp>
 #include <opm/input/eclipse/EclipseState/Tables/Regdims.hpp>
 #include <opm/input/eclipse/EclipseState/Tables/Tabdims.hpp>
 
@@ -635,6 +636,7 @@ class Runspec {
 
     const Tracers& tracers() const;
     const Geochem& geochem() const;
+    const PorosityModel& porosityModel() const noexcept;
     bool compositionalMode() const;
     std::size_t numComps() const;
     std::size_t maxGasPlantTables() const;
@@ -683,6 +685,7 @@ class Runspec {
         serializer(m_mechsolver);
         serializer(m_biof);
         serializer(m_geochem);
+        serializer(m_porosity_model);
     }
 
 private:
@@ -703,6 +706,7 @@ class Runspec {
     MechSolver m_mechsolver{};
     Tracers m_tracers{};
     Geochem m_geochem{};
+    PorosityModel m_porosity_model{};
     std::size_t m_comps = 0;
     std::size_t m_max_gas_plant_tables = 0;
     bool m_co2storage{false};

opm/input/eclipse/share/keywords/000_Eclipse100/N/NODPPM

@@ -1,7 +1,6 @@
 {
   "name": "NODPPM",
   "sections": [
-    "RUNSPEC",
-    "GRID"
+    "RUNSPEC"
   ]
 }

opm/input/eclipse/share/keywords/000_Eclipse100/S/SIGMAV

@@ -4,6 +4,7 @@
     "GRID"
   ],
   "data": {
-    "value_type": "DOUBLE"
+    "value_type": "DOUBLE",
+    "dimension": "1/Length*Length"
   }
 }

tests/parser/FieldPropsTests.cpp

@@ -155,6 +155,51 @@ PERMX
 
 
 
+BOOST_AUTO_TEST_CASE(SigmaVFieldProps) {
+    // SIGMAV is a per-cell array (JSON schema: "data", no "size") -- fits FieldProps'
+    // GRID::double_keywords registry exactly like PORO/PERMX.
+    std::string deck_string_sigmav = R"(
+GRID
+
+PORO
+   8*0.10 /
+
+SIGMAV
+  8*0.12 /
+)";
+    EclipseGrid grid(EclipseGrid(2,2,2));
+    Deck deck_sigmav = Parser{}.parseString(deck_string_sigmav);
+    FieldPropsManager fpm_sigmav(deck_sigmav, Phases{true, true, false}, grid, TableManager());
+
+    BOOST_CHECK(fpm_sigmav.has_double("SIGMAV"));
+    const auto& sigmav = fpm_sigmav.get_double("SIGMAV");
+    BOOST_CHECK_EQUAL(sigmav.size(), grid.getNumActive());
+    for (const auto& value : sigmav) {
+        BOOST_CHECK_CLOSE(value, 0.12, 1e-10);
+    }
+}
+
+BOOST_AUTO_TEST_CASE(SigmaDeckScalar) {
+    // SIGMA is a single global scalar (JSON schema: "size": 1), not a per-cell array like
+    // SIGMAV (schema: "data", no "size"). FieldProps::GRID::double_keywords is strictly a
+    // per-cell array registry: verify_deck_data() in FieldProps.cpp requires
+    // deck_data.size() == box.size() * num_value unconditionally, with no broadcast path for
+    // a single supplied value. SIGMA is intentionally not registered there; in dual-continuum
+    // runs applyDualPorosityScalars broadcasts its value into the SIGMAV carrier instead.
+    // This test proves the raw parse of a size-1 keyword works standalone.
+    std::string deck_string_sigma = R"(
+GRID
+
+SIGMA
+  0.12 /
+)";
+    auto deck_sigma = Parser{}.parseString(deck_string_sigma);
+    BOOST_CHECK(deck_sigma.hasKeyword("SIGMA"));
+    const auto sigma_value = deck_sigma["SIGMA"].back().getRecord(0).getItem(0).get<double>(0);
+    BOOST_CHECK_CLOSE(sigma_value, 0.12, 1e-10);
+}
+
+
 BOOST_AUTO_TEST_CASE(CreateFieldPropsForActnum) {
     std::string deck_string = R"(
 GRID