From 3f6c4922b0c16044bba31df4e05de9045da5b2c2 Mon Sep 17 00:00:00 2001 From: Arne Morten Kvarving Date: Tue, 8 Nov 2022 07:09:51 +0100 Subject: [PATCH] move updatePrimaryVariablesNewton to StandardWellPrimaryVariables --- opm/simulators/wells/StandardWellEval.cpp | 51 ------------------- opm/simulators/wells/StandardWellEval.hpp | 4 -- .../wells/StandardWellPrimaryVariables.cpp | 49 ++++++++++++++++++ .../wells/StandardWellPrimaryVariables.hpp | 15 ++++-- opm/simulators/wells/StandardWell_impl.hpp | 2 +- 5 files changed, 60 insertions(+), 61 deletions(-) diff --git a/opm/simulators/wells/StandardWellEval.cpp b/opm/simulators/wells/StandardWellEval.cpp index 4ee80ccb8..3b7fb6b96 100644 --- a/opm/simulators/wells/StandardWellEval.cpp +++ b/opm/simulators/wells/StandardWellEval.cpp @@ -107,57 +107,6 @@ computeAccumWell() } } -template -void -StandardWellEval:: -updatePrimaryVariablesNewton(const BVectorWell& dwells, - [[maybe_unused]] const double dFLimit, - const double dBHPLimit) const -{ - const std::vector old_primary_variables = primary_variables_.value_; - - // for injectors, very typical one of the fractions will be one, and it is easy to get zero value - // fractions. not sure what is the best way to handle it yet, so we just use 1.0 here - [[maybe_unused]] const double relaxation_factor_fractions = - (baseif_.isProducer()) ? this->primary_variables_.relaxationFactorFractionsProducer(old_primary_variables, dwells) - : 1.0; - - // update the second and third well variable (The flux fractions) - - if constexpr (has_wfrac_variable) { - const int sign2 = dwells[0][WFrac] > 0 ? 1: -1; - const double dx2_limited = sign2 * std::min(std::abs(dwells[0][WFrac] * relaxation_factor_fractions), dFLimit); - // primary_variables_[WFrac] = old_primary_variables[WFrac] - dx2_limited; - primary_variables_.value_[WFrac] = old_primary_variables[WFrac] - dx2_limited; - } - - if constexpr (has_gfrac_variable) { - const int sign3 = dwells[0][GFrac] > 0 ? 1: -1; - const double dx3_limited = sign3 * std::min(std::abs(dwells[0][GFrac] * relaxation_factor_fractions), dFLimit); - primary_variables_.value_[GFrac] = old_primary_variables[GFrac] - dx3_limited; - } - - if constexpr (Indices::enableSolvent) { - const int sign4 = dwells[0][SFrac] > 0 ? 1: -1; - const double dx4_limited = sign4 * std::min(std::abs(dwells[0][SFrac]) * relaxation_factor_fractions, dFLimit); - primary_variables_.value_[SFrac] = old_primary_variables[SFrac] - dx4_limited; - } - - this->primary_variables_.processFractions(); - - // updating the total rates Q_t - const double relaxation_factor_rate = this->relaxationFactorRate(old_primary_variables, dwells[0][WQTotal]); - primary_variables_.value_[WQTotal] = old_primary_variables[WQTotal] - dwells[0][WQTotal] * relaxation_factor_rate; - - // updating the bottom hole pressure - { - const int sign1 = dwells[0][Bhp] > 0 ? 1: -1; - const double dx1_limited = sign1 * std::min(std::abs(dwells[0][Bhp]), std::abs(old_primary_variables[Bhp]) * dBHPLimit); - // 1e5 to make sure bhp will not be below 1bar - primary_variables_.value_[Bhp] = std::max(old_primary_variables[Bhp] - dx1_limited, 1e5); - } -} - template ConvergenceReport StandardWellEval:: diff --git a/opm/simulators/wells/StandardWellEval.hpp b/opm/simulators/wells/StandardWellEval.hpp index 66940bd0b..71a16e017 100644 --- a/opm/simulators/wells/StandardWellEval.hpp +++ b/opm/simulators/wells/StandardWellEval.hpp @@ -118,10 +118,6 @@ protected: void updateWellStateFromPrimaryVariables(WellState& well_state, DeferredLogger& deferred_logger) const; - void updatePrimaryVariablesNewton(const BVectorWell& dwells, - const double dFLimit, - const double dBHPLimit) const; - mutable PrimaryVariables primary_variables_; //!< Primary variables for well // the saturations in the well bore under surface conditions at the beginning of the time step diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.cpp b/opm/simulators/wells/StandardWellPrimaryVariables.cpp index d6f30e560..58bd2acd7 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.cpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.cpp @@ -215,6 +215,55 @@ updatePolyMW(const BVectorWell& dwells) } } +template +void StandardWellPrimaryVariables:: +updateNewton(const BVectorWell& dwells, + [[maybe_unused]] const double dFLimit, + const double dBHPLimit) +{ + const std::vector old_primary_variables = value_; + + // for injectors, very typical one of the fractions will be one, and it is easy to get zero value + // fractions. not sure what is the best way to handle it yet, so we just use 1.0 here + [[maybe_unused]] const double relaxation_factor_fractions = + well_.isProducer() ? relaxationFactorFractionsProducer(old_primary_variables, dwells) + : 1.0; + + // update the second and third well variable (The flux fractions) + + if constexpr (has_wfrac_variable) { + const int sign2 = dwells[0][WFrac] > 0 ? 1: -1; + const double dx2_limited = sign2 * std::min(std::abs(dwells[0][WFrac] * relaxation_factor_fractions), dFLimit); + // primary_variables_[WFrac] = old_primary_variables[WFrac] - dx2_limited; + value_[WFrac] = old_primary_variables[WFrac] - dx2_limited; + } + + if constexpr (has_gfrac_variable) { + const int sign3 = dwells[0][GFrac] > 0 ? 1: -1; + const double dx3_limited = sign3 * std::min(std::abs(dwells[0][GFrac] * relaxation_factor_fractions), dFLimit); + value_[GFrac] = old_primary_variables[GFrac] - dx3_limited; + } + + if constexpr (Indices::enableSolvent) { + const int sign4 = dwells[0][SFrac] > 0 ? 1: -1; + const double dx4_limited = sign4 * std::min(std::abs(dwells[0][SFrac]) * relaxation_factor_fractions, dFLimit); + value_[SFrac] = old_primary_variables[SFrac] - dx4_limited; + } + + this->processFractions(); + + // updating the total rates Q_t + const double relaxation_factor_rate = StandardWellGeneric::relaxationFactorRate(old_primary_variables, dwells[0][WQTotal]); + value_[WQTotal] = old_primary_variables[WQTotal] - dwells[0][WQTotal] * relaxation_factor_rate; + + // updating the bottom hole pressure + const int sign1 = dwells[0][Bhp] > 0 ? 1: -1; + const double dx1_limited = sign1 * std::min(std::abs(dwells[0][Bhp]), + std::abs(old_primary_variables[Bhp]) * dBHPLimit); + // 1e5 to make sure bhp will not be below 1bar + value_[Bhp] = std::max(old_primary_variables[Bhp] - dx1_limited, 1e5); +} + template void StandardWellPrimaryVariables:: copyToWellState(WellState& well_state, diff --git a/opm/simulators/wells/StandardWellPrimaryVariables.hpp b/opm/simulators/wells/StandardWellPrimaryVariables.hpp index 9b0a9d22e..f9ffe1549 100644 --- a/opm/simulators/wells/StandardWellPrimaryVariables.hpp +++ b/opm/simulators/wells/StandardWellPrimaryVariables.hpp @@ -111,7 +111,12 @@ public: //! \brief Copy values from well state. void update(const WellState& well_state, DeferredLogger& deferred_logger); - //! \brief Update polymer molecular weight values from solution vector. + //! \brief Update values from newton update vector. + void updateNewton(const BVectorWell& dwells, + const double dFLimit, + const double dBHPLimit); + + //! \brief Update polymer molecular weight values from newton update vector. void updatePolyMW(const BVectorWell& dwells); //! \brief Copy values to well state. @@ -129,18 +134,18 @@ public: //! \brief Returns scaled rate for a component. EvalWell getQs(const int compIdx) const; - //! \brief Handle non-reasonable fractions due to numerical overshoot. - void processFractions(); - +private: //! \brief Calculate a relaxation factor for producers. //! \details To avoid overshoot of the fractions which might result in negative rates. double relaxationFactorFractionsProducer(const std::vector& primary_variables, const BVectorWell& dwells) const; -private: //! \brief Returns volume fraction for a component. EvalWell volumeFraction(const unsigned compIdx) const; + //! \brief Handle non-reasonable fractions due to numerical overshoot. + void processFractions(); + const WellInterfaceIndices& well_; //!< Reference to well interface //! \brief Total number of the well equations and primary variables. diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index a723b1560..e75f606cc 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -957,7 +957,7 @@ namespace Opm { const double dFLimit = this->param_.dwell_fraction_max_; const double dBHPLimit = this->param_.dbhp_max_rel_; - this->StdWellEval::updatePrimaryVariablesNewton(dwells, dFLimit, dBHPLimit); + this->primary_variables_.updateNewton(dwells, dFLimit, dBHPLimit); updateExtraPrimaryVariables(dwells);