diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 8f20f07d1..aaa414c8d 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -106,4 +106,5 @@ list (APPEND PUBLIC_HEADER_FILES opm/autodiff/SimulatorIncompTwophaseAd.hpp opm/autodiff/TransportSolverTwophaseAd.hpp opm/autodiff/WellDensitySegmented.hpp + opm/autodiff/WellStateFullyImplicitBlackoil.hpp ) diff --git a/examples/sim_fibo_ad.cpp b/examples/sim_fibo_ad.cpp index cc305a78d..74353a503 100644 --- a/examples/sim_fibo_ad.cpp +++ b/examples/sim_fibo_ad.cpp @@ -39,7 +39,7 @@ #include #include -#include +#include #include #include @@ -205,7 +205,7 @@ try std::cout << "\n\n================ Starting main simulation loop ===============\n" << std::flush; - WellState well_state; + WellStateFullyImplicitBlackoil well_state; Opm::TimeMapPtr timeMap(new Opm::TimeMap(newParserDeck)); SimulatorTimer simtimer; @@ -248,7 +248,7 @@ try SimulatorReport rep; // With a deck, we may have more epochs etc. - WellState well_state; + WellStateFullyImplicitBlackoil well_state; int step = 0; SimulatorTimer simtimer; // Use timer for last epoch to obtain total time. diff --git a/examples/test_implicit_ad.cpp b/examples/test_implicit_ad.cpp index 0c1e78580..11558b396 100644 --- a/examples/test_implicit_ad.cpp +++ b/examples/test_implicit_ad.cpp @@ -23,6 +23,7 @@ #include #include #include +#include #include #include @@ -37,7 +38,6 @@ #include #include -#include #include #include @@ -110,7 +110,7 @@ try initStateBasic(*g, props0, param, 0.0, state); initBlackoilSurfvol(*g, props0, state); - Opm::WellState well_state; + Opm::WellStateFullyImplicitBlackoil well_state; well_state.init(wells.get(), state); solver.step(1.0, state, well_state); diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.cpp b/opm/autodiff/FullyImplicitBlackoilSolver.cpp index f3f686712..b0b53c44e 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.cpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.cpp @@ -24,12 +24,12 @@ #include #include #include +#include #include #include #include #include -#include #include #include @@ -223,7 +223,7 @@ namespace { FullyImplicitBlackoilSolver:: step(const double dt, BlackoilState& x , - WellState& xw) + WellStateFullyImplicitBlackoil& xw) { const V pvdt = geo_.poreVolume() / dt; @@ -331,7 +331,7 @@ namespace { FullyImplicitBlackoilSolver::SolutionState FullyImplicitBlackoilSolver::constantState(const BlackoilState& x, - const WellState& xw) + const WellStateFullyImplicitBlackoil& xw) { const int nc = grid_.number_of_cells; const int np = x.numPhases(); @@ -427,7 +427,7 @@ namespace { FullyImplicitBlackoilSolver::SolutionState FullyImplicitBlackoilSolver::variableState(const BlackoilState& x, - const WellState& xw) + const WellStateFullyImplicitBlackoil& xw) { const int nc = grid_.number_of_cells; const int np = x.numPhases(); @@ -616,7 +616,7 @@ namespace { FullyImplicitBlackoilSolver:: assemble(const V& pvdt, const BlackoilState& x , - const WellState& xw ) + const WellStateFullyImplicitBlackoil& xw ) { // Create the primary variables. const SolutionState state = variableState(x, xw); @@ -878,7 +878,7 @@ namespace { void FullyImplicitBlackoilSolver::updateState(const V& dx, BlackoilState& state, - WellState& well_state) + WellStateFullyImplicitBlackoil& well_state) { const int np = fluid_.numPhases(); const int nc = grid_.number_of_cells; diff --git a/opm/autodiff/FullyImplicitBlackoilSolver.hpp b/opm/autodiff/FullyImplicitBlackoilSolver.hpp index 65d888026..985a3a571 100644 --- a/opm/autodiff/FullyImplicitBlackoilSolver.hpp +++ b/opm/autodiff/FullyImplicitBlackoilSolver.hpp @@ -33,7 +33,7 @@ namespace Opm { class RockCompressibility; class LinearSolverInterface; class BlackoilState; - class WellState; + class WellStateFullyImplicitBlackoil; /// A fully implicit solver for the black-oil problem. @@ -76,7 +76,7 @@ namespace Opm { void step(const double dt , BlackoilState& state , - WellState& wstate); + WellStateFullyImplicitBlackoil& wstate); private: // Types and enums @@ -148,11 +148,11 @@ namespace Opm { // Private methods. SolutionState constantState(const BlackoilState& x, - const WellState& xw); + const WellStateFullyImplicitBlackoil& xw); SolutionState variableState(const BlackoilState& x, - const WellState& xw); + const WellStateFullyImplicitBlackoil& xw); void computeAccum(const SolutionState& state, @@ -160,14 +160,14 @@ namespace Opm { void assemble(const V& dtpv, - const BlackoilState& x , - const WellState& xw ); + const BlackoilState& x, + const WellStateFullyImplicitBlackoil& xw); V solveJacobianSystem() const; void updateState(const V& dx, BlackoilState& state, - WellState& well_state); + WellStateFullyImplicitBlackoil& well_state); std::vector computePressures(const SolutionState& state) const; diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoil.cpp b/opm/autodiff/SimulatorFullyImplicitBlackoil.cpp index 00246f90f..b827cabce 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoil.cpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoil.cpp @@ -29,6 +29,7 @@ #include #include #include +#include #include #include @@ -48,7 +49,6 @@ #include #include -#include #include #include @@ -77,7 +77,7 @@ namespace Opm SimulatorReport run(SimulatorTimer& timer, BlackoilState& state, - WellState& well_state); + WellStateFullyImplicitBlackoil& well_state); private: // Data. @@ -127,7 +127,7 @@ namespace Opm SimulatorReport SimulatorFullyImplicitBlackoil::run(SimulatorTimer& timer, BlackoilState& state, - WellState& well_state) + WellStateFullyImplicitBlackoil& well_state) { return pimpl_->run(timer, state, well_state); } @@ -198,7 +198,7 @@ namespace Opm std::copy(d.begin(), d.end(), std::ostream_iterator(file, "\n")); } } - static void outputWellStateMatlab(const Opm::WellState& well_state, + static void outputWellStateMatlab(const Opm::WellStateFullyImplicitBlackoil& well_state, const int step, const std::string& output_dir) { @@ -312,10 +312,10 @@ namespace Opm SimulatorReport SimulatorFullyImplicitBlackoil::Impl::run(SimulatorTimer& timer, BlackoilState& state, - WellState& well_state) + WellStateFullyImplicitBlackoil& well_state) { - eclipseWriter_.writeInit(timer, state, well_state); - eclipseWriter_.writeTimeStep(timer, state, well_state); + eclipseWriter_.writeInit(timer, state, well_state.basicWellState()); + eclipseWriter_.writeTimeStep(timer, state, well_state.basicWellState()); // Initialisation. std::vector porevol; @@ -384,13 +384,9 @@ namespace Opm // Optionally, check if well controls are satisfied. if (check_well_controls_) { - Opm::computePhaseFlowRatesPerWell(*wells_, - well_state.perfRates(), - fractional_flows, - well_resflows_phase); std::cout << "Checking well conditions." << std::endl; // For testing we set surface := reservoir - well_control_passed = wells_manager_.conditionsMet(well_state.bhp(), well_resflows_phase, well_resflows_phase); + well_control_passed = wells_manager_.conditionsMet(well_state.bhp(), well_state.wellRates(), well_state.wellRates()); ++well_control_iteration; if (!well_control_passed && well_control_iteration > max_well_control_iterations_) { OPM_THROW(std::runtime_error, "Could not satisfy well conditions in " << max_well_control_iterations_ << " tries."); @@ -426,7 +422,7 @@ namespace Opm // write an output file for later inspection if (output_) { - eclipseWriter_.writeTimeStep(timer, state, well_state); + eclipseWriter_.writeTimeStep(timer, state, well_state.basicWellState()); } // advance to next timestep before reporting at this location diff --git a/opm/autodiff/SimulatorFullyImplicitBlackoil.hpp b/opm/autodiff/SimulatorFullyImplicitBlackoil.hpp index 2bb62de41..26b91b4e3 100644 --- a/opm/autodiff/SimulatorFullyImplicitBlackoil.hpp +++ b/opm/autodiff/SimulatorFullyImplicitBlackoil.hpp @@ -37,7 +37,7 @@ namespace Opm class LinearSolverInterface; class SimulatorTimer; class BlackoilState; - class WellState; + class WellStateFullyImplicitBlackoil; struct SimulatorReport; /// Class collecting all necessary components for a two-phase simulation. @@ -84,7 +84,7 @@ namespace Opm /// \return simulation report, with timing data SimulatorReport run(SimulatorTimer& timer, BlackoilState& state, - WellState& well_state); + WellStateFullyImplicitBlackoil& well_state); private: class Impl; diff --git a/opm/autodiff/WellStateFullyImplicitBlackoil.hpp b/opm/autodiff/WellStateFullyImplicitBlackoil.hpp new file mode 100644 index 000000000..f2430b1b8 --- /dev/null +++ b/opm/autodiff/WellStateFullyImplicitBlackoil.hpp @@ -0,0 +1,101 @@ +/* + Copyright 2014 SINTEF ICT, Applied Mathematics. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED +#define OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED + + +#include +#include +#include +#include +#include + +namespace Opm +{ + + /// The state of a set of wells, tailored for use by the fully + /// implicit blackoil simulator. + class WellStateFullyImplicitBlackoil + { + public: + /// Allocate and initialize if wells is non-null. Also tries + /// to give useful initial values to the bhp(), wellRates() + /// and perfPhaseRates() fields, depending on controls + template + void init(const Wells* wells, const State& state) + { + if (wells == 0) { + return; + } + + // We use the WellState::init() function to do bhp and well rates init. + // The alternative would be to copy that function wholesale. + basic_well_state_.init(wells, state); + + // Initialize perfphaserates_, which must be done here. + const int nw = wells->number_of_wells; + const int np = wells->number_of_phases; + const int nperf = wells->well_connpos[nw]; + perfphaserates_.resize(nperf * np, 0.0); + for (int w = 0; w < nw; ++w) { + assert((wells->type[w] == INJECTOR) || (wells->type[w] == PRODUCER)); + const WellControls* ctrl = wells->ctrls[w]; + if (well_controls_well_is_shut(ctrl)) { + // Shut well: perfphaserates_ are all zero. + } else { + // Open well: Initialize perfphaserates_ to well + // rates divided by the number of perforations. + const int num_perf_this_well = wells->well_connpos[w + 1] - wells->well_connpos[w]; + for (int perf = wells->well_connpos[w]; perf < wells->well_connpos[w + 1]; ++perf) { + for (int p = 0; p < np; ++p) { + perfphaserates_[np*perf + p] = wellRates()[np*w + p] / double(num_perf_this_well); + } + } + } + } + } + + /// One bhp pressure per well. + std::vector& bhp() { return basic_well_state_.bhp(); } + const std::vector& bhp() const { return basic_well_state_.bhp(); } + + /// One rate per well and phase. + std::vector& wellRates() { return basic_well_state_.wellRates(); } + const std::vector& wellRates() const { return basic_well_state_.wellRates(); } + + /// One rate per phase and well connection. + std::vector& perfPhaseRates() { return perfphaserates_; } + const std::vector& perfPhaseRates() const { return perfphaserates_; } + + /// For interfacing with functions that take a WellState. + const WellState& basicWellState() const + { + return basic_well_state_; + } + + private: + WellState basic_well_state_; + std::vector perfphaserates_; + }; + +} // namespace Opm + + +#endif // OPM_WELLSTATEFULLYIMPLICITBLACKOIL_HEADER_INCLUDED