diff --git a/examples/wells_example.cpp b/examples/wells_example.cpp index 68b54ce5..2cc9bd8e 100644 --- a/examples/wells_example.cpp +++ b/examples/wells_example.cpp @@ -3,8 +3,11 @@ #include #include #include - +#include +#include +#include #include "opm/core/newwells.h" +#include "opm/core/grid.h" int main(int argc, char** argv) { @@ -12,7 +15,7 @@ int main(int argc, char** argv) { using namespace Opm; ParameterGroup parameters( argc, argv, false ); std::string file_name = parameters.getDefault("inputdeck", "data.data"); - + // Read input file EclipseGridParser parser(file_name); std::cout << "Done!" << std::endl; @@ -22,6 +25,18 @@ int main(int argc, char** argv) { // Finally handle the wells WellsManager wells(parser, *grid.c_grid(), NULL); + std::vector global_cells(grid.c_grid()->global_cell, grid.c_grid()->global_cell + grid.c_grid()->number_of_cells); + + std::cout << "ahoi" << std::endl; + double gravity[3] = {0.0, 0.0, parameters.getDefault("gravity", 0.0)}; + IncompPropertiesFromDeck incomp_properties(parser, global_cells); + std::cout << "there" << std::endl; + + LinearSolverUmfpack umfpack_solver; + std::cout << "here" << std::endl; + IncompTpfa pressure_solver(*grid.c_grid(), incomp_properties.permeability(), + gravity, umfpack_solver, wells.c_wells()); + return 0; } diff --git a/opm/core/WellCollection.cpp b/opm/core/WellCollection.cpp index e2568f55..9d1b9f27 100644 --- a/opm/core/WellCollection.cpp +++ b/opm/core/WellCollection.cpp @@ -41,11 +41,11 @@ namespace Opm } std::tr1::shared_ptr child; - for (int i = 0; i < roots_.size(); ++i) { + for (size_t i = 0; i < roots_.size(); ++i) { if (roots_[i]->name() == child_name) { child = roots_[i]; // We've found a new parent to the previously thought root, need to remove it - for(int j = i; j < roots_.size() - 1; ++j) { + for(size_t j = i; j < roots_.size() - 1; ++j) { roots_[j] = roots_[j+1]; } @@ -63,13 +63,19 @@ namespace Opm } parent_as_group->addChild(child); - + if(child->isLeafNode()) { + leaf_nodes_.push_back(child); + } + } + + const std::vector >& WellCollection::getLeafNodes() const { + return leaf_nodes_; } WellsGroupInterface* WellCollection::findNode(std::string name) { - for (int i = 0; i < roots_.size(); i++) { + for (size_t i = 0; i < roots_.size(); i++) { WellsGroupInterface* result = roots_[i]->findGroup(name); if (result) { return result; diff --git a/opm/core/WellCollection.hpp b/opm/core/WellCollection.hpp index 53ef008e..a4cd782c 100644 --- a/opm/core/WellCollection.hpp +++ b/opm/core/WellCollection.hpp @@ -38,9 +38,14 @@ namespace Opm void addChild(std::string child, std::string parent, const EclipseGridParser& deck); + const std::vector >& getLeafNodes() const; private: // To account for the possibility of a forest std::vector > roots_; + + // This will be used to traverse the bottom nodes. + std::vector > leaf_nodes_; + WellsGroupInterface* findNode(std::string name); }; diff --git a/opm/core/WellsGroup.cpp b/opm/core/WellsGroup.cpp index 4b36acf8..17fa8a4d 100644 --- a/opm/core/WellsGroup.cpp +++ b/opm/core/WellsGroup.cpp @@ -31,13 +31,17 @@ namespace Opm : WellsGroupInterface(name, prod_spec, inj_spec) { } + + bool WellsGroupInterface::isLeafNode() const { + return false; + } WellsGroupInterface* WellsGroup::findGroup(std::string name_of_node) { if (name() == name_of_node) { return this; } else { - for (int i = 0; i < children_.size(); i++) { + for (size_t i = 0; i < children_.size(); i++) { WellsGroupInterface* result = children_[i]->findGroup(name_of_node); if (result) { return result; @@ -68,6 +72,10 @@ namespace Opm } } + + bool WellNode::isLeafNode() const { + return true; + } surface_component toSurfaceComponent(std::string type) { @@ -171,7 +179,7 @@ namespace Opm bool isWell = false; if (deck.hasField("WELSPECS")) { WELSPECS wspecs = deck.getWELSPECS(); - for (int i = 0; i < wspecs.welspecs.size(); i++) { + for (size_t i = 0; i < wspecs.welspecs.size(); i++) { if (wspecs.welspecs[i].name_ == name) { isWell = true; break; @@ -185,7 +193,7 @@ namespace Opm InjectionSpecification injection_specification; if (deck.hasField("WCONINJE")) { WCONINJE wconinje = deck.getWCONINJE(); - for (int i = 0; i < wconinje.wconinje.size(); i++) { + for (size_t i = 0; i < wconinje.wconinje.size(); i++) { if (wconinje.wconinje[i].well_ == name) { WconinjeLine line = wconinje.wconinje[i]; injection_specification.BHP_limit_ = line.BHP_limit_; @@ -199,7 +207,7 @@ namespace Opm ProductionSpecification production_specification; if (deck.hasField("WCONPROD")) { WCONPROD wconprod = deck.getWCONPROD(); - for (int i = 0; i < wconprod.wconprod.size(); i++) { + for (size_t i = 0; i < wconprod.wconprod.size(); i++) { if (wconprod.wconprod[i].well_ == name) { WconprodLine line = wconprod.wconprod[i]; production_specification.BHP_limit_ = line.BHP_limit_; @@ -216,7 +224,7 @@ namespace Opm InjectionSpecification injection_specification; if (deck.hasField("GCONINJE")) { GCONINJE gconinje = deck.getGCONINJE(); - for (int i = 0; i < gconinje.gconinje.size(); i++) { + for (size_t i = 0; i < gconinje.gconinje.size(); i++) { if (gconinje.gconinje[i].group_ == name) { GconinjeLine line = gconinje.gconinje[i]; injection_specification.injector_type_ = toSurfaceComponent(line.injector_type_); @@ -229,7 +237,7 @@ namespace Opm ProductionSpecification production_specification; if (deck.hasField("GCONPROD")) { GCONPROD gconprod = deck.getGCONPROD(); - for (int i = 0; i < gconprod.gconprod.size(); i++) { + for (size_t i = 0; i < gconprod.gconprod.size(); i++) { if (gconprod.gconprod[i].group_ == name) { GconprodLine line = gconprod.gconprod[i]; production_specification.oil_max_rate_ = line.oil_max_rate_; diff --git a/opm/core/WellsGroup.hpp b/opm/core/WellsGroup.hpp index 9d5131ce..7a9ad14d 100644 --- a/opm/core/WellsGroup.hpp +++ b/opm/core/WellsGroup.hpp @@ -17,9 +17,17 @@ namespace Opm InjectionSpecification inj_spec); virtual ~WellsGroupInterface(); + /// The unique identifier for the well or well group. const std::string& name(); + + /// Production specifications for the well or well group. const ProductionSpecification& prodSpec() const; + + /// Injection specifications for the well or well group. const InjectionSpecification& injSpec() const; + + /// \returns true if the object is a leaf node (WellNode), false otherwise. + virtual bool isLeafNode() const; /// \returns the pointer to the WellsGroupInterface with the given name. NULL if /// the name is not found.a @@ -51,6 +59,8 @@ namespace Opm InjectionSpecification inj_spec); virtual WellsGroupInterface* findGroup(std::string name_of_node); + + virtual bool isLeafNode() const; }; diff --git a/opm/core/WellsManager.cpp b/opm/core/WellsManager.cpp index 36bc12b5..2d7fb5b2 100644 --- a/opm/core/WellsManager.cpp +++ b/opm/core/WellsManager.cpp @@ -505,22 +505,22 @@ namespace Opm } } - WellCollection well_collection; + WellCollection wells_; if (deck.hasField("GRUPTREE")) { std::cout << "Found gruptree" << std::endl; const GRUPTREE& gruptree = deck.getGRUPTREE(); std::map::const_iterator it = gruptree.tree.begin(); for( ; it != gruptree.tree.end(); ++it) { - well_collection.addChild(it->first, it->second, deck); + wells_.addChild(it->first, it->second, deck); } } if(deck.hasField("WELSPECS")) { WELSPECS welspecs = deck.getWELSPECS(); - for(int i = 0; i < welspecs.welspecs.size(); ++i) { + for(size_t i = 0; i < welspecs.welspecs.size(); ++i) { WelspecsLine line = welspecs.welspecs[i]; - well_collection.addChild(line.name_, line.group_, deck); + wells_.addChild(line.name_, line.group_, deck); } } diff --git a/opm/core/WellsManager.hpp b/opm/core/WellsManager.hpp index e17eaeb7..36070dd6 100644 --- a/opm/core/WellsManager.hpp +++ b/opm/core/WellsManager.hpp @@ -19,7 +19,7 @@ #ifndef OPM_WELLSMANAGER_HEADER_INCLUDED #define OPM_WELLSMANAGER_HEADER_INCLUDED - +#include struct Wells; struct UnstructuredGrid; @@ -61,6 +61,8 @@ namespace Opm WellsManager& operator=(const WellsManager& other); // The managed Wells. Wells* w_; + + WellCollection well_collection_; }; diff --git a/opm/core/pressure/IncompTpfa.cpp b/opm/core/pressure/IncompTpfa.cpp index 9417b7eb..8ec414d2 100644 --- a/opm/core/pressure/IncompTpfa.cpp +++ b/opm/core/pressure/IncompTpfa.cpp @@ -38,14 +38,18 @@ namespace Opm /// and N == g.number_of_cells. /// \param[in] gravity Gravity vector. If nonzero, the array should /// have D elements. + /// \param[in] wells The wells argument. Will be used in solution, + /// is ignored if NULL IncompTpfa::IncompTpfa(const UnstructuredGrid& g, const double* permeability, const double* gravity, - const LinearSolverInterface& linsolver) + const LinearSolverInterface& linsolver, + const struct Wells* wells) : grid_(g), linsolver_(linsolver), htrans_(g.cell_facepos[ g.number_of_cells ]), - trans_ (g.number_of_faces) + trans_ (g.number_of_faces), + wells_(wells) { UnstructuredGrid* gg = const_cast(&grid_); tpfa_htrans_compute(gg, permeability, &htrans_[0]); @@ -110,7 +114,7 @@ namespace Opm } } - ifs_tpfa_forces F = { NULL, NULL, NULL, NULL, NULL }; + ifs_tpfa_forces F = { NULL, NULL, wells_, NULL, NULL }; if (! src.empty()) { F.src = &src[0]; } F.bc = bcs; @@ -172,7 +176,7 @@ namespace Opm } } - ifs_tpfa_forces F = { NULL, NULL, NULL, NULL, NULL }; + ifs_tpfa_forces F = { NULL, NULL, wells_, NULL, NULL }; if (! src.empty()) { F.src = &src[0]; } F.bc = bcs; diff --git a/opm/core/pressure/IncompTpfa.hpp b/opm/core/pressure/IncompTpfa.hpp index 58b53383..01eb6627 100644 --- a/opm/core/pressure/IncompTpfa.hpp +++ b/opm/core/pressure/IncompTpfa.hpp @@ -25,6 +25,7 @@ struct UnstructuredGrid; struct ifs_tpfa_data; +struct Wells; struct FlowBoundaryConditions; namespace Opm @@ -47,10 +48,13 @@ namespace Opm /// and N == g.number_of_cells. /// \param[in] gravity Gravity vector. If nonzero, the array should /// have D elements. + /// \param[in] wells The wells argument. Will be used in solution, + /// is ignored if NULL IncompTpfa(const UnstructuredGrid& g, const double* permeability, const double* gravity, - const LinearSolverInterface& linsolver); + const LinearSolverInterface& linsolver, + const struct Wells* wells = 0); /// Destructor. ~IncompTpfa(); @@ -113,7 +117,8 @@ namespace Opm ::std::vector trans_ ; ::std::vector gpress_; ::std::vector gpress_omegaweighted_; - + + const struct Wells* wells_; struct ifs_tpfa_data* h_; }; diff --git a/opm/core/utility/writeVtkData.cpp b/opm/core/utility/writeVtkData.cpp index 63c8284c..e331e3e5 100644 --- a/opm/core/utility/writeVtkData.cpp +++ b/opm/core/utility/writeVtkData.cpp @@ -22,6 +22,7 @@ #include #include #include +#include #include #include #include