From e22234d6a60a7556f94f3d85d532c1f73ebe884f Mon Sep 17 00:00:00 2001 From: plgbrts Date: Thu, 3 Jul 2025 15:50:18 +0200 Subject: [PATCH] water viscosity calculation for brine + thermal --- .../blackoilpvt/WaterPvtThermal.cpp | 42 +++++++++++++------ .../blackoilpvt/WaterPvtThermal.hpp | 30 +++++++++---- 2 files changed, 52 insertions(+), 20 deletions(-) diff --git a/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.cpp b/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.cpp index 9b207c0453f..402297ce758 100644 --- a/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.cpp +++ b/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.cpp @@ -56,6 +56,7 @@ initFromState(const EclipseState& eclState, const Schedule& schedule) enableJouleThomson_ = tables.WatJT().size() > 0; enableThermalViscosity_ = tables.hasTables("WATVISCT"); enableInternalEnergy_ = tables.hasTables("SPECHEAT"); + enableBrineViscosity_ = eclState.runspec().phases().active(Phase::BRINE); unsigned regions = isothermalPvt_->numRegions(); setNumRegions(regions); @@ -115,13 +116,30 @@ initFromState(const EclipseState& eclState, const Schedule& schedule) const auto& watvisctTables = tables.getWatvisctTables(); const auto& viscrefTables = tables.getViscrefTable(); - const auto& pvtwTables = tables.getPvtwTable(); - - if (pvtwTables.size() != regions) { - OPM_THROW(std::runtime_error, - fmt::format("Table sizes mismatch. PVTW: {}, numRegions: {}\n", - pvtwTables.size(), regions)); + if (enableBrineViscosity_) { + const auto& pvtwSaltTables = tables.getPvtwSaltTables(); + if (pvtwSaltTables.size() != regions) { + OPM_THROW(std::runtime_error, + fmt::format("Table sizes mismatch. PVTWSALT: {}, numRegions: {}\n", + pvtwSaltTables.size(), regions)); + } + for (unsigned regionIdx = 0; regionIdx < regions; ++ regionIdx) { + referenceSaltConcentration_[regionIdx] = pvtwSaltTables[regionIdx].getReferenceSaltConcentrationValue(); + } } + else { + const auto& pvtwTables = tables.getPvtwTable(); + if (pvtwTables.size() != regions) { + OPM_THROW(std::runtime_error, + fmt::format("Table sizes mismatch. PVTW: {}, numRegions: {}\n", + pvtwTables.size(), regions)); + } + for (unsigned regionIdx = 0; regionIdx < regions; ++ regionIdx) { + pvtwViscosity_[regionIdx] = pvtwTables[regionIdx].viscosity; + pvtwViscosibility_[regionIdx] = pvtwTables[regionIdx].viscosibility; + } + } + if (watvisctTables.size() != regions) { OPM_THROW(std::runtime_error, fmt::format("Table sizes mismatch. WATVISCT: {}, numRegions: {}\n", @@ -140,11 +158,6 @@ initFromState(const EclipseState& eclState, const Schedule& schedule) viscrefPress_[regionIdx] = viscrefTables[regionIdx].reference_pressure; } - - for (unsigned regionIdx = 0; regionIdx < regions; ++ regionIdx) { - pvtwViscosity_[regionIdx] = pvtwTables[regionIdx].viscosity; - pvtwViscosibility_[regionIdx] = pvtwTables[regionIdx].viscosibility; - } } if (enableInternalEnergy_) { @@ -189,6 +202,7 @@ setNumRegions(std::size_t numRegions) pvtwCompressibility_.resize(numRegions); pvtwViscosity_.resize(numRegions); pvtwViscosibility_.resize(numRegions); + referenceSaltConcentration_.resize(numRegions); viscrefPress_.resize(numRegions); watvisctCurves_.resize(numRegions); watdentRefTemp_.resize(numRegions); @@ -222,12 +236,14 @@ operator==(const WaterPvtThermal& data) const this->pvtwCompressibility() == data.pvtwCompressibility() && this->pvtwViscosity() == data.pvtwViscosity() && this->pvtwViscosibility() == data.pvtwViscosibility() && + this->referenceSaltConcentration() == data.referenceSaltConcentration() && this->watvisctCurves() == data.watvisctCurves() && this->internalEnergyCurves() == data.internalEnergyCurves() && this->enableThermalDensity() == data.enableThermalDensity() && this->enableJouleThomson() == data.enableJouleThomson() && this->enableThermalViscosity() == data.enableThermalViscosity() && - this->enableInternalEnergy() == data.enableInternalEnergy(); + this->enableInternalEnergy() == data.enableInternalEnergy() && + this->enableBrineViscosity() == data.enableBrineViscosity(); } template @@ -252,12 +268,14 @@ operator=(const WaterPvtThermal& data) pvtwCompressibility_ = data.pvtwCompressibility_; pvtwViscosity_ = data.pvtwViscosity_; pvtwViscosibility_ = data.pvtwViscosibility_; + referenceSaltConcentration_ = data.referenceSaltConcentration_; watvisctCurves_ = data.watvisctCurves_; internalEnergyCurves_ = data.internalEnergyCurves_; enableThermalDensity_ = data.enableThermalDensity_; enableJouleThomson_ = data.enableJouleThomson_; enableThermalViscosity_ = data.enableThermalViscosity_; enableInternalEnergy_ = data.enableInternalEnergy_; + enableBrineViscosity_ = data.enableBrineViscosity_; return *this; } diff --git a/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.hpp b/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.hpp index 55b89481a11..97fac6a8ac4 100644 --- a/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.hpp +++ b/opm/material/fluidsystems/blackoilpvt/WaterPvtThermal.hpp @@ -68,12 +68,14 @@ class WaterPvtThermal const std::vector& pvtwCompressibility, const std::vector& pvtwViscosity, const std::vector& pvtwViscosibility, + const std::vector& referenceSaltConcentration, const std::vector& watvisctCurves, const std::vector& internalEnergyCurves, bool enableThermalDensity, bool enableJouleThomson, bool enableThermalViscosity, - bool enableInternalEnergy) + bool enableInternalEnergy, + bool enableBrineViscosity) : isothermalPvt_(isothermalPvt) , viscrefPress_(viscrefPress) , watdentRefTemp_(watdentRefTemp) @@ -86,12 +88,14 @@ class WaterPvtThermal , pvtwCompressibility_(pvtwCompressibility) , pvtwViscosity_(pvtwViscosity) , pvtwViscosibility_(pvtwViscosibility) + , referenceSaltConcentration_(referenceSaltConcentration) , watvisctCurves_(watvisctCurves) , internalEnergyCurves_(internalEnergyCurves) , enableThermalDensity_(enableThermalDensity) , enableJouleThomson_(enableJouleThomson) , enableThermalViscosity_(enableThermalViscosity) , enableInternalEnergy_(enableInternalEnergy) + , enableBrineViscosity_(enableBrineViscosity) { } WaterPvtThermal(const WaterPvtThermal& data) @@ -227,14 +231,16 @@ class WaterPvtThermal if (!enableThermalViscosity()) { return isothermalMu; } - - Scalar x = -pvtwViscosibility_[regionIdx] * (viscrefPress_[regionIdx] - - pvtwRefPress_[regionIdx]); - Scalar muRef = pvtwViscosity_[regionIdx] / (1.0 + x + 0.5 * x * x); - - // compute the viscosity deviation due to temperature + const auto& muWatvisct = watvisctCurves_[regionIdx].eval(temperature, true); - return isothermalMu * muWatvisct / muRef; + if (enableBrineViscosity()) { + auto muRef = isothermalPvt_->viscosity(regionIdx, temperature, Evaluation(viscrefPress_[regionIdx]), Rsw, Evaluation(referenceSaltConcentration_[regionIdx])); + return isothermalMu * muWatvisct / muRef; + } else { + Scalar x = -pvtwViscosibility_[regionIdx] * (viscrefPress_[regionIdx] - pvtwRefPress_[regionIdx]); + Scalar muRef = pvtwViscosity_[regionIdx] / (1.0 + x + 0.5 * x * x); + return isothermalMu * muWatvisct / muRef; + } } /*! @@ -404,6 +410,9 @@ class WaterPvtThermal const std::vector& pvtwViscosibility() const { return pvtwViscosibility_; } + const std::vector& referenceSaltConcentration() const + { return referenceSaltConcentration_; } + const std::vector& watvisctCurves() const { return watvisctCurves_; } @@ -413,6 +422,9 @@ class WaterPvtThermal bool enableInternalEnergy() const { return enableInternalEnergy_; } + bool enableBrineViscosity() const + { return enableBrineViscosity_; } + const std::vector& watJTRefPres() const { return watJTRefPres_; } @@ -443,6 +455,7 @@ class WaterPvtThermal std::vector pvtwCompressibility_{}; std::vector pvtwViscosity_{}; std::vector pvtwViscosibility_{}; + std::vector referenceSaltConcentration_{}; std::vector watvisctCurves_{}; @@ -454,6 +467,7 @@ class WaterPvtThermal bool enableJouleThomson_{false}; bool enableThermalViscosity_{false}; bool enableInternalEnergy_{false}; + bool enableBrineViscosity_{false}; }; } // namespace Opm