From 22dadf5928c397ff45316ce719acbda86f38979d Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Tue, 4 Feb 2025 09:43:24 +0100 Subject: [PATCH] component indices are now int inactive components are represented using negative values in the index(traits) class so we should not be using unsigned values in the fluidsystem while at it make activeToCanonicalComponentIndex usable at compile time --- opm/models/blackoil/blackoilindices.hh | 4 ++-- opm/models/blackoil/blackoilonephaseindices.hh | 8 +++++--- opm/models/blackoil/blackoiltwophaseindices.hh | 8 ++++---- .../wells/MultisegmentWellPrimaryVariables.cpp | 4 ++-- .../wells/MultisegmentWellPrimaryVariables.hpp | 4 ++-- opm/simulators/wells/RatioCalculator.cpp | 6 +++--- opm/simulators/wells/RatioCalculator.hpp | 12 ++++++------ opm/simulators/wells/StandardWellConnections.cpp | 2 +- .../wells/StandardWellPrimaryVariables.cpp | 4 ++-- .../wells/StandardWellPrimaryVariables.hpp | 2 +- opm/simulators/wells/WellInterfaceIndices.cpp | 2 +- opm/simulators/wells/WellInterfaceIndices.hpp | 2 +- 12 files changed, 30 insertions(+), 28 deletions(-) diff --git a/opm/models/blackoil/blackoilindices.hh b/opm/models/blackoil/blackoilindices.hh index 8350373d0..53f4ed6bb 100644 --- a/opm/models/blackoil/blackoilindices.hh +++ b/opm/models/blackoil/blackoilindices.hh @@ -94,10 +94,10 @@ struct BlackOilIndices numEnergy + numFoam + numBrine + numMICPs; //! \brief returns the index of "active" component - static constexpr unsigned canonicalToActiveComponentIndex(unsigned compIdx) + static constexpr int canonicalToActiveComponentIndex(const int compIdx) { return compIdx; } - static constexpr unsigned activeToCanonicalComponentIndex(unsigned compIdx) + static constexpr int activeToCanonicalComponentIndex(const int compIdx) { return compIdx; } //////// diff --git a/opm/models/blackoil/blackoilonephaseindices.hh b/opm/models/blackoil/blackoilonephaseindices.hh index c2abedb02..863c8b00e 100644 --- a/opm/models/blackoil/blackoilonephaseindices.hh +++ b/opm/models/blackoil/blackoilonephaseindices.hh @@ -30,6 +30,8 @@ #include +#include + namespace Opm { /*! @@ -179,15 +181,15 @@ struct BlackOilOnePhaseIndices ////////////////////// //! \brief returns the index of "active" component - static constexpr unsigned canonicalToActiveComponentIndex(unsigned /*compIdx*/) + static constexpr int canonicalToActiveComponentIndex(const int /*compIdx*/) { return 0; } - static unsigned activeToCanonicalComponentIndex([[maybe_unused]] unsigned compIdx) + static constexpr int activeToCanonicalComponentIndex([[maybe_unused]] const int compIdx) { // assumes canonical oil = 0, water = 1, gas = 2; - assert(compIdx == 0); + constexpr_assert(compIdx == 0); if (gasEnabled) { return 2; } else if (waterEnabled) { diff --git a/opm/models/blackoil/blackoiltwophaseindices.hh b/opm/models/blackoil/blackoiltwophaseindices.hh index e4eefa50d..9d112a148 100644 --- a/opm/models/blackoil/blackoiltwophaseindices.hh +++ b/opm/models/blackoil/blackoiltwophaseindices.hh @@ -181,7 +181,7 @@ struct BlackOilTwoPhaseIndices ////////////////////// //! \brief returns the index of "active" component - static constexpr unsigned canonicalToActiveComponentIndex(const unsigned compIdx) + static constexpr int canonicalToActiveComponentIndex(const int compIdx) { // assumes canonical oil = 0, water = 1, gas = 2; if (!gasEnabled) { @@ -201,10 +201,10 @@ struct BlackOilTwoPhaseIndices return compIdx - 1; } - static unsigned activeToCanonicalComponentIndex(unsigned compIdx) + static constexpr int activeToCanonicalComponentIndex(const int compIdx) { // assumes canonical oil = 0, water = 1, gas = 2; - assert(compIdx < 2); + constexpr_assert(compIdx < 2); if (!gasEnabled) { // oil = 0, water = 1 return compIdx; @@ -212,7 +212,7 @@ struct BlackOilTwoPhaseIndices // oil = 0, gas = 1 return compIdx * 2; } else { - assert(!oilEnabled); + constexpr_assert(!oilEnabled); } // water = 0, gas = 1; diff --git a/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp b/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp index 1c2a99507..f3d410b24 100644 --- a/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp +++ b/opm/simulators/wells/MultisegmentWellPrimaryVariables.cpp @@ -515,7 +515,7 @@ template typename MultisegmentWellPrimaryVariables::EvalWell MultisegmentWellPrimaryVariables:: volumeFraction(const int seg, - const unsigned compIdx) const + const int compIdx) const { if (has_wfrac_variable && compIdx == Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)) { return evaluation_[seg][WFrac]; @@ -584,7 +584,7 @@ typename MultisegmentWellPrimaryVariables::EvalWell MultisegmentWellPrimaryVariables:: getSegmentRateUpwinding(const int seg, const int seg_upwind, - const std::size_t comp_idx) const + const int comp_idx) const { // the result will contain the derivative with respect to WQTotal in segment seg, // and the derivatives with respect to WFrac GFrac in segment seg_upwind. diff --git a/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp b/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp index 3dd05a61c..951d34b4c 100644 --- a/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp +++ b/opm/simulators/wells/MultisegmentWellPrimaryVariables.hpp @@ -115,7 +115,7 @@ public: //! \brief Returns upwinding rate for a component in a segment. EvalWell getSegmentRateUpwinding(const int seg, const int seg_upwind, - const std::size_t comp_idx) const; + const int comp_idx) const; //! \brief Get bottomhole pressure. EvalWell getBhp() const; @@ -154,7 +154,7 @@ private: //! \brief Returns volume fraction for component in a segment. EvalWell volumeFraction(const int seg, - const unsigned compIdx) const; + const int compIdx) const; //! \brief The values for the primary variables //! \details Based on different solution strategies, the wells can have different primary variables diff --git a/opm/simulators/wells/RatioCalculator.cpp b/opm/simulators/wells/RatioCalculator.cpp index d64990b7c..72d8ea5cd 100644 --- a/opm/simulators/wells/RatioCalculator.cpp +++ b/opm/simulators/wells/RatioCalculator.cpp @@ -51,9 +51,9 @@ namespace Opm { template RatioCalculator:: -RatioCalculator(unsigned gasCompIdx, - unsigned oilCompIdx, - unsigned waterCompIdx, +RatioCalculator(int gasCompIdx, + int oilCompIdx, + int waterCompIdx, std::string_view name) : gasComp_{gasCompIdx} , oilComp_(oilCompIdx) diff --git a/opm/simulators/wells/RatioCalculator.hpp b/opm/simulators/wells/RatioCalculator.hpp index 8121c3bb1..985881b88 100644 --- a/opm/simulators/wells/RatioCalculator.hpp +++ b/opm/simulators/wells/RatioCalculator.hpp @@ -39,9 +39,9 @@ class RatioCalculator public: using Scalar = decltype(getValue(Value{})); - RatioCalculator(unsigned gasCompIdx, - unsigned oilCompIdx, - unsigned waterCompIdx, + RatioCalculator(int gasCompIdx, + int oilCompIdx, + int waterCompIdx, std::string_view name); void disOilVapWatVolumeRatio(Value& volumeRatio, @@ -91,9 +91,9 @@ public: const bool isProducer) const; private: - unsigned gasComp_; - unsigned oilComp_; - unsigned waterComp_; + int gasComp_; + int oilComp_; + int waterComp_; std::string name_; }; diff --git a/opm/simulators/wells/StandardWellConnections.cpp b/opm/simulators/wells/StandardWellConnections.cpp index 3b4c17f20..d333ee0fc 100644 --- a/opm/simulators/wells/StandardWellConnections.cpp +++ b/opm/simulators/wells/StandardWellConnections.cpp @@ -731,7 +731,7 @@ connectionRateFoam(const std::vector& cq_s, } case Phase::SOLVENT: { if constexpr (Indices::enableSolvent) - return static_cast(Indices::contiSolventEqIdx); + return Indices::contiSolventEqIdx; else OPM_DEFLOG_THROW(std::runtime_error, "Foam transport phase is SOLVENT but SOLVENT is not activated.", deferred_logger); } diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.cpp b/opm/simulators/wells/StandardWellPrimaryVariables.cpp index 17410e2f3..3f7369ef9 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.cpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.cpp @@ -446,7 +446,7 @@ copyToWellStatePolyMW(WellState& well_state) const template typename StandardWellPrimaryVariables::EvalWell StandardWellPrimaryVariables:: -volumeFraction(const unsigned compIdx) const +volumeFraction(const int compIdx) const { if (FluidSystem::numActivePhases() == 1) { return EvalWell(numWellEq_ + Indices::numEq, 1.0); @@ -456,7 +456,7 @@ volumeFraction(const unsigned compIdx) const return evaluation_[GFrac]; } - if (Indices::enableSolvent && compIdx == (unsigned)Indices::contiSolventEqIdx) { + if (Indices::enableSolvent && compIdx == Indices::contiSolventEqIdx) { return evaluation_[SFrac]; } diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.hpp b/opm/simulators/wells/StandardWellPrimaryVariables.hpp index f482af12e..7c0c551ac 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.hpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.hpp @@ -157,7 +157,7 @@ private: DeferredLogger& deferred_logger) const; //! \brief Returns volume fraction for a component. - EvalWell volumeFraction(const unsigned compIdx) const; + EvalWell volumeFraction(const int compIdx) const; //! \brief Handle non-reasonable fractions due to numerical overshoot. void processFractions(); diff --git a/opm/simulators/wells/WellInterfaceIndices.cpp b/opm/simulators/wells/WellInterfaceIndices.cpp index 1150bdc18..983318488 100644 --- a/opm/simulators/wells/WellInterfaceIndices.cpp +++ b/opm/simulators/wells/WellInterfaceIndices.cpp @@ -78,7 +78,7 @@ flowPhaseToModelCompIdx(const int phaseIdx) const template int WellInterfaceIndices:: -modelCompIdxToFlowCompIdx(const unsigned compIdx) const +modelCompIdxToFlowCompIdx(const int compIdx) const { const auto& pu = this->phaseUsage(); if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx) == compIdx) diff --git a/opm/simulators/wells/WellInterfaceIndices.hpp b/opm/simulators/wells/WellInterfaceIndices.hpp index 6ac6138e1..16146c9a9 100644 --- a/opm/simulators/wells/WellInterfaceIndices.hpp +++ b/opm/simulators/wells/WellInterfaceIndices.hpp @@ -41,7 +41,7 @@ public: using ModelParameters = typename WellInterfaceFluidSystem::ModelParameters; int flowPhaseToModelCompIdx(const int phaseIdx) const; - int modelCompIdxToFlowCompIdx(const unsigned compIdx) const; + int modelCompIdxToFlowCompIdx(const int compIdx) const; Scalar scalingFactor(const int phaseIdx) const; template