diff --git a/opm/autodiff/StandardWellsDense.hpp b/opm/autodiff/StandardWellsDense.hpp index 33e9a550b..f42e7e8a7 100644 --- a/opm/autodiff/StandardWellsDense.hpp +++ b/opm/autodiff/StandardWellsDense.hpp @@ -363,6 +363,8 @@ enum WellVariablePositions { const int well_index, WellState& xw) const; + bool wellHasTHPConstraints(const int well_index) const; + }; diff --git a/opm/autodiff/StandardWellsDense_impl.hpp b/opm/autodiff/StandardWellsDense_impl.hpp index cd47fb334..beff786e6 100644 --- a/opm/autodiff/StandardWellsDense_impl.hpp +++ b/opm/autodiff/StandardWellsDense_impl.hpp @@ -1739,6 +1739,11 @@ namespace Opm { well_potentials.resize(nw * np, 0.0); for (int w = 0; w < nw; ++w) { + + bool is_thp_determined = wellHasTHPConstraints(w); + + if (!is_thp_determined) { + // bhp needs to be determined for the well potential calculation // There can be more than one BHP/THP constraints. // TODO: there is an option to ignore the THP limit when calculating well potentials, @@ -2751,4 +2756,23 @@ namespace Opm { } + + + + template + bool + StandardWellsDense:: + wellHasTHPConstraints(const int well_index) const + { + const WellType& well_type = wells().type[well_index]; + const WellControls* well_control = wells().ctrls[well_index]; + const int nwc = well_controls_get_num(well_control); + for (int ctrl_index = 0; ctrl_index < nwc; ++ctrl_index) { + if (well_controls_iget_type(well_control, ctrl_index) == THP) { + return true; + } + } + return false; + } + } // namespace Opm