diff --git a/examples/spu_2p.cpp b/examples/spu_2p.cpp index a6128d2b..467438f1 100644 --- a/examples/spu_2p.cpp +++ b/examples/spu_2p.cpp @@ -510,6 +510,7 @@ main(int argc, char** argv) } else { computeTotalMobility(*props, allcells, state.saturation(), totmob); } + std::vector empty_vector_for_wells; pressure_timer.start(); if (rock_comp->isActive()) { rc.resize(num_cells); @@ -521,8 +522,10 @@ main(int argc, char** argv) rc[cell] = rock_comp->rockComp(state.pressure()[cell]); } state.pressure() = initial_pressure; + + psolver.solve(totmob, omega, src, bcs.c_bcs(), porevol, rc, simtimer.currentStepLength(), - state.pressure(), state.faceflux()); + state.pressure(), state.faceflux(), empty_vector_for_wells, empty_vector_for_wells); double max_change = 0.0; for (int cell = 0; cell < num_cells; ++cell) { max_change = std::max(max_change, std::fabs(state.pressure()[cell] - prev_pressure[cell])); @@ -534,7 +537,8 @@ main(int argc, char** argv) } computePorevolume(*grid->c_grid(), *props, *rock_comp, state.pressure(), porevol); } else { - psolver.solve(totmob, omega, src, bcs.c_bcs(), state.pressure(), state.faceflux()); + psolver.solve(totmob, omega, src, bcs.c_bcs(), state.pressure(), state.faceflux(), empty_vector_for_wells, + empty_vector_for_wells); } pressure_timer.stop(); double pt = pressure_timer.secsSinceStart(); diff --git a/opm/core/ProductionSpecification.cpp b/opm/core/ProductionSpecification.cpp index abc7d061..1c23b2ff 100644 --- a/opm/core/ProductionSpecification.cpp +++ b/opm/core/ProductionSpecification.cpp @@ -10,7 +10,8 @@ namespace Opm oil_max_rate_(-1.0), water_production_target_(-1.0), fluid_volume_max_rate_(-1.0), - BHP_limit_(-1.0) + BHP_limit_(-1.0), + guide_rate_(-1.0) { } diff --git a/opm/core/ProductionSpecification.hpp b/opm/core/ProductionSpecification.hpp index 16243126..afd20ac0 100644 --- a/opm/core/ProductionSpecification.hpp +++ b/opm/core/ProductionSpecification.hpp @@ -18,6 +18,11 @@ namespace Opm { NONE_P, RATE, WELL }; + + enum GuideRateType + { + OIL, RAT + }; ProductionSpecification(); @@ -30,6 +35,8 @@ namespace Opm double BHP_limit_; double gas_max_rate_; double liquid_max_rate_; + double guide_rate_; + GuideRateType guide_rate_type_; }; } diff --git a/opm/core/WellsGroup.cpp b/opm/core/WellsGroup.cpp index 99ed02e1..e7968f21 100644 --- a/opm/core/WellsGroup.cpp +++ b/opm/core/WellsGroup.cpp @@ -97,9 +97,7 @@ namespace Opm } } - - - + return true; } WellsGroupInterface* WellNode::findGroup(std::string name_of_node) diff --git a/opm/core/WellsManager.cpp b/opm/core/WellsManager.cpp index 03c60fae..46100e24 100644 --- a/opm/core/WellsManager.cpp +++ b/opm/core/WellsManager.cpp @@ -527,6 +527,11 @@ namespace Opm // We know that getLeafNodes() is ordered the same way as they're indexed in w_ node->setWellsPointer(w_, i); } + + // Set the guide rates: + if(deck.hasField("WCONPROD")) { + WCONPROD wconprod = deck.getWCONPROD(); + } } diff --git a/opm/core/pressure/IncompTpfa.cpp b/opm/core/pressure/IncompTpfa.cpp index 8ec414d2..540c1ebd 100644 --- a/opm/core/pressure/IncompTpfa.cpp +++ b/opm/core/pressure/IncompTpfa.cpp @@ -91,12 +91,18 @@ namespace Opm /// If null, noflow conditions are assumed. /// \param[out] pressure Will contain N cell-pressure values. /// \param[out] faceflux Will contain F signed face flux values. + /// \param[out] well_bhp Will contain bhp values for each well passed + /// in the constructor. + /// \param[out] well_rate Will contain rate values for each well passed + /// in the constructor. void IncompTpfa::solve(const std::vector& totmob, const std::vector& omega, const std::vector& src, const FlowBoundaryConditions* bcs, std::vector& pressure, - std::vector& faceflux) + std::vector& faceflux, + std::vector& well_bhp, + std::vector& well_rate) { UnstructuredGrid* gg = const_cast(&grid_); tpfa_eff_trans_compute(gg, &totmob[0], &htrans_[0], &trans_[0]); @@ -150,6 +156,10 @@ namespace Opm /// \param[in] dt Timestep. /// \param[out] pressure Will contain N cell-pressure values. /// \param[out] faceflux Will contain F signed face flux values. + /// \param[out] well_bhp Will contain bhp values for each well passed + /// in the constructor + /// \param[out] well_rate Will contain rate values for each well passed + /// in the constructor void IncompTpfa::solve(const std::vector& totmob, const std::vector& omega, const std::vector& src, @@ -158,7 +168,9 @@ namespace Opm const std::vector& rock_comp, const double dt, std::vector& pressure, - std::vector& faceflux) + std::vector& faceflux, + std::vector& well_bhp, + std::vector& well_rate) { UnstructuredGrid* gg = const_cast(&grid_); tpfa_eff_trans_compute(gg, &totmob[0], &htrans_[0], &trans_[0]); diff --git a/opm/core/pressure/IncompTpfa.hpp b/opm/core/pressure/IncompTpfa.hpp index 01eb6627..23105ce1 100644 --- a/opm/core/pressure/IncompTpfa.hpp +++ b/opm/core/pressure/IncompTpfa.hpp @@ -72,12 +72,18 @@ namespace Opm /// If null, noflow conditions are assumed. /// \param[out] pressure Will contain N cell-pressure values. /// \param[out] faceflux Will contain F signed face flux values. + /// \param[out] well_bhp Will contain bhp values for each well passed + /// in the constructor + /// \param[out] well_rate Will contain rate values for each well passed + /// in the constructor void solve(const std::vector& totmob, const std::vector& omega, const std::vector& src, const FlowBoundaryConditions* bcs, std::vector& pressure, - std::vector& faceflux); + std::vector& faceflux, + std::vector& well_bhp, + std::vector& well_rate); /// Assemble and solve pressure system with rock compressibility (assumed constant per cell). /// \param[in] totmob Must contain N total mobility values (one per cell). @@ -97,6 +103,10 @@ namespace Opm /// \param[in] dt Timestep. /// \param[out] pressure Will contain N cell-pressure values. /// \param[out] faceflux Will contain F signed face flux values. + /// \param[out] well_bhp Will contain bhp values for each well passed + /// in the constructor + /// \param[out] well_rate Will contain rate values for each well passed + /// in the constructor void solve(const std::vector& totmob, const std::vector& omega, const std::vector& src, @@ -105,7 +115,9 @@ namespace Opm const std::vector& rock_comp, const double dt, std::vector& pressure, - std::vector& faceflux); + std::vector& faceflux, + std::vector& well_bhp, + std::vector& well_rate); /// Expose read-only reference to internal half-transmissibility. const ::std::vector& getHalfTrans() const { return htrans_; }