From c5a0ea75247f5cf350348eacd5766f4ae90b4e55 Mon Sep 17 00:00:00 2001 From: Andreas Lauser Date: Thu, 26 Jan 2017 17:36:00 +0100 Subject: [PATCH] do not explicitly pass the permeability to the well model anymore this information is already part of the EclipseState. The reason why this should IMO be avoided is that this enforces an implementation (ordering of the permeability matrices) the simulator on the well model. If this needs to be done for performance reasons, IMO it would be smarter to pass an array of matrices, instead of passing a raw array of doubles. I doubt that this is necessary, though: completing the full Norne deck takes about 0.25 seconds longer on my machine, that's substantially less than 0.1% of the total runtime. in order to avoid code duplication, the permeability extraction function of the RockFromDeck class is now made a public static function and used as an implementation detail of the WellsManager. finally, the permfield_valid_ attribute is removed from the RockFromDeck class because this data was unused and not accessible via the class' public API. --- examples/compute_tof.cpp | 2 +- examples/wells_example.cpp | 2 +- opm/core/props/rock/RockFromDeck.cpp | 32 +++++++++++++------------- opm/core/props/rock/RockFromDeck.hpp | 22 ++++++++++++++---- opm/core/wells/WellsManager.cpp | 5 ++-- opm/core/wells/WellsManager.hpp | 5 +--- opm/core/wells/WellsManager_impl.hpp | 16 +++++++++---- tests/test_stoppedwells.cpp | 4 ++-- tests/test_wellsmanager.cpp | 18 +++++++-------- tests/wells_manager_data.data | 2 ++ tests/wells_manager_data_wellSTOP.data | 3 +++ tests/wells_stopped.data | 2 ++ 12 files changed, 68 insertions(+), 45 deletions(-) diff --git a/examples/compute_tof.cpp b/examples/compute_tof.cpp index 89a5013d..63a74651 100644 --- a/examples/compute_tof.cpp +++ b/examples/compute_tof.cpp @@ -182,7 +182,7 @@ try // Rock and fluid init IncompPropertiesSinglePhase props(deck, eclipseState, grid); // Wells init. - WellsManager wells_manager(eclipseState , 0, grid, props.permeability()); + WellsManager wells_manager(eclipseState , 0, grid); std::shared_ptr my_wells(clone_wells(wells_manager.c_wells()), destroy_wells); setBhpWells(*my_wells); diff --git a/examples/wells_example.cpp b/examples/wells_example.cpp index af0c2031..48f68c8b 100644 --- a/examples/wells_example.cpp +++ b/examples/wells_example.cpp @@ -49,7 +49,7 @@ try RockCompressibility rock_comp(eclipseState); // Finally handle the wells - WellsManager wells(eclipseState , 0 , *grid.c_grid(), incomp_properties.permeability()); + WellsManager wells(eclipseState , 0 , *grid.c_grid()); double gravity[3] = {0.0, 0.0, parameters.getDefault("gravity", 0.0)}; Opm::LinearSolverFactory linsolver(parameters); diff --git a/opm/core/props/rock/RockFromDeck.cpp b/opm/core/props/rock/RockFromDeck.cpp index 20eec193..84d9c876 100644 --- a/opm/core/props/rock/RockFromDeck.cpp +++ b/opm/core/props/rock/RockFromDeck.cpp @@ -74,8 +74,7 @@ namespace Opm RockFromDeck::RockFromDeck(std::size_t number_of_cells) : porosity_(number_of_cells, 0), // full permeability tensor in 3D stores 9 scalars - permeability_(number_of_cells*9, 0.0), - permfield_valid_(number_of_cells, false) + permeability_(number_of_cells*9, 0.0) { } @@ -84,10 +83,13 @@ namespace Opm const int* cart_dims) { assignPorosity(eclState, number_of_cells, global_cell); - permfield_valid_.assign(number_of_cells, false); const double perm_threshold = 0.0; // Maybe turn into parameter? - assignPermeability(eclState, number_of_cells, global_cell, cart_dims, - perm_threshold); + extractInterleavedPermeability(eclState, + number_of_cells, + global_cell, + cart_dims, + perm_threshold, + permeability_); } void RockFromDeck::assignPorosity(const Opm::EclipseState& eclState, @@ -105,11 +107,12 @@ namespace Opm } } - void RockFromDeck::assignPermeability(const Opm::EclipseState& eclState, - int number_of_cells, - const int* global_cell, - const int* cartdims, - double perm_threshold) + void RockFromDeck::extractInterleavedPermeability(const Opm::EclipseState& eclState, + const int number_of_cells, + const int* global_cell, + const int* cartdims, + const double perm_threshold, + std::vector& permeability) { const int dim = 3; const int nc = number_of_cells; @@ -117,7 +120,7 @@ namespace Opm assert(cartdims[0]*cartdims[1]*cartdims[2] > 0); static_cast(cartdims); // Squash warning in release mode. - permeability_.assign(dim * dim * nc, 0.0); + permeability.assign(dim * dim * nc, 0.0); std::vector tensor; tensor.reserve(6); @@ -132,7 +135,6 @@ namespace Opm assert (! tensor.empty()); { int off = 0; - for (int c = 0; c < nc; ++c, off += dim*dim) { // SharedPermTensor K(dim, dim, &permeability_[off]); int kix = 0; @@ -146,16 +148,14 @@ namespace Opm // values in the resulting array are the same // in either order when viewed contiguously // because fillTensor() enforces symmetry. - permeability_[off + (i + dim*j)] = + permeability[off + (i + dim*j)] = tensor[kmap[kix]][c]; } // K(i,i) = std::max(K(i,i), perm_threshold); - double& kii = permeability_[off + i*(dim + 1)]; + double& kii = permeability[off + i*(dim + 1)]; kii = std::max(kii, perm_threshold); } - - permfield_valid_[c] = std::vector::value_type(1); } } } diff --git a/opm/core/props/rock/RockFromDeck.hpp b/opm/core/props/rock/RockFromDeck.hpp index 00e74958..c96fa1dc 100644 --- a/opm/core/props/rock/RockFromDeck.hpp +++ b/opm/core/props/rock/RockFromDeck.hpp @@ -77,15 +77,27 @@ namespace Opm return &permeability_[0]; } + /// Convert the permeabilites for the logically Cartesian grid in EclipseState to + /// an array of size number_of_cells*dim*dim for the compressed array. + /// \param eclState The EclipseState (processed deck) produced by the opm-parser code + /// \param number_of_cells The number of cells in the grid. + /// \param global_cell The mapping fom local to global cell indices. + /// global_cell[i] is the corresponding global index of i. + /// \param cart_dims The size of the underlying cartesian grid. + /// \param perm_threshold The threshold for permeability + /// \param permeability The result array + static + void extractInterleavedPermeability(const Opm::EclipseState& eclState, + const int number_of_cells, + const int* global_cell, + const int* cart_dims, + const double perm_threshold, + std::vector& permeability); + private: void assignPorosity(const Opm::EclipseState& eclState, int number_of_cells, const int* global_cell); - void assignPermeability(const Opm::EclipseState& eclState, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - double perm_threshold); std::vector porosity_; std::vector permeability_; diff --git a/opm/core/wells/WellsManager.cpp b/opm/core/wells/WellsManager.cpp index 7e139652..79397216 100644 --- a/opm/core/wells/WellsManager.cpp +++ b/opm/core/wells/WellsManager.cpp @@ -328,8 +328,7 @@ namespace Opm /// Construct wells from deck. WellsManager::WellsManager(const Opm::EclipseState& eclipseState, const size_t timeStep, - const UnstructuredGrid& grid, - const double* permeability) + const UnstructuredGrid& grid) : w_(0), is_parallel_run_(false) { std::vector dummy_well_potentials; @@ -340,7 +339,7 @@ namespace Opm UgGridHelpers::globalCell(grid), UgGridHelpers::cartDims(grid), UgGridHelpers::dimensions(grid), UgGridHelpers::cell2Faces(grid), UgGridHelpers::beginFaceCentroids(grid), - permeability, dummy_list_econ_limited, dummy_well_potentials, + dummy_list_econ_limited, dummy_well_potentials, std::unordered_set()); } diff --git a/opm/core/wells/WellsManager.hpp b/opm/core/wells/WellsManager.hpp index 3593202a..adcbef84 100644 --- a/opm/core/wells/WellsManager.hpp +++ b/opm/core/wells/WellsManager.hpp @@ -91,7 +91,6 @@ namespace Opm int dimensions, const F2C& f2c, FC begin_face_centroids, - const double* permeability, const DynamicListEconLimited& list_econ_limited, bool is_parallel_run=false, const std::vector& well_potentials=std::vector(), @@ -99,8 +98,7 @@ namespace Opm WellsManager(const Opm::EclipseState& eclipseState, const size_t timeStep, - const UnstructuredGrid& grid, - const double* permeability); + const UnstructuredGrid& grid); /// Destructor. ~WellsManager(); @@ -163,7 +161,6 @@ namespace Opm int dimensions, const C2F& cell_to_faces, FC begin_face_centroids, - const double* permeability, const DynamicListEconLimited& list_econ_limited, const std::vector& well_potentials, const std::unordered_set& deactivated_wells); diff --git a/opm/core/wells/WellsManager_impl.hpp b/opm/core/wells/WellsManager_impl.hpp index 30324bac..cac0df18 100644 --- a/opm/core/wells/WellsManager_impl.hpp +++ b/opm/core/wells/WellsManager_impl.hpp @@ -4,6 +4,7 @@ #include #include #include +#include #include #include #include @@ -311,7 +312,6 @@ WellsManager(const Opm::EclipseState& eclipseState, int dimensions, const C2F& cell_to_faces, FC begin_face_centroids, - const double* permeability, const DynamicListEconLimited& list_econ_limited, bool is_parallel_run, const std::vector& well_potentials, @@ -320,7 +320,7 @@ WellsManager(const Opm::EclipseState& eclipseState, { init(eclipseState, timeStep, number_of_cells, global_cell, cart_dims, dimensions, - cell_to_faces, begin_face_centroids, permeability, list_econ_limited, well_potentials, deactivated_wells); + cell_to_faces, begin_face_centroids, list_econ_limited, well_potentials, deactivated_wells); } /// Construct wells from deck. @@ -334,7 +334,6 @@ WellsManager::init(const Opm::EclipseState& eclipseState, int dimensions, const C2F& cell_to_faces, FC begin_face_centroids, - const double* permeability, const DynamicListEconLimited& list_econ_limited, const std::vector& well_potentials, const std::unordered_set& deactivated_wells) @@ -392,13 +391,22 @@ WellsManager::init(const Opm::EclipseState& eclipseState, } } + + std::vector interleavedPerm; + RockFromDeck::extractInterleavedPermeability(eclipseState, + number_of_cells, + global_cell, + cart_dims, + 0.0, + interleavedPerm); + createWellsFromSpecs(wells, timeStep, cell_to_faces, cart_dims, begin_face_centroids, dimensions, dz, well_names, well_data, well_names_to_index, - pu, cartesian_to_compressed, permeability, ntg, + pu, cartesian_to_compressed, interleavedPerm.data(), ntg, wells_on_proc, deactivated_wells, list_econ_limited); setupWellControls(wells, timeStep, well_names, pu, wells_on_proc, list_econ_limited); diff --git a/tests/test_stoppedwells.cpp b/tests/test_stoppedwells.cpp index c1ea918d..73d847c0 100644 --- a/tests/test_stoppedwells.cpp +++ b/tests/test_stoppedwells.cpp @@ -58,7 +58,7 @@ BOOST_AUTO_TEST_CASE(TestStoppedWells) // Both wells are open in the first schedule step { - Opm::WellsManager wellsManager(eclipseState , 0 , *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState , 0 , *gridManager.c_grid()); const Wells* wells = wellsManager.c_wells(); const struct WellControls* ctrls0 = wells->ctrls[0]; const struct WellControls* ctrls1 = wells->ctrls[1]; @@ -78,7 +78,7 @@ BOOST_AUTO_TEST_CASE(TestStoppedWells) // The injector is stopped { - Opm::WellsManager wellsManager(eclipseState , 1 , *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState , 1 , *gridManager.c_grid()); const Wells* wells = wellsManager.c_wells(); const struct WellControls* ctrls0 = wells->ctrls[0]; const struct WellControls* ctrls1 = wells->ctrls[1]; diff --git a/tests/test_wellsmanager.cpp b/tests/test_wellsmanager.cpp index 774a9146..38d6c68e 100644 --- a/tests/test_wellsmanager.cpp +++ b/tests/test_wellsmanager.cpp @@ -185,19 +185,19 @@ BOOST_AUTO_TEST_CASE(New_Constructor_Works) { Opm::GridManager gridManager(eclipseState.getInputGrid()); { - Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid()); wells_static_check(wellsManager.c_wells()); check_controls_epoch0(wellsManager.c_wells()->ctrls); } { - Opm::WellsManager wellsManager(eclipseState, 1, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState, 1, *gridManager.c_grid()); wells_static_check(wellsManager.c_wells()); check_controls_epoch1(wellsManager.c_wells()->ctrls); } { - Opm::WellsManager wellsManager(eclipseState, 3, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState, 3, *gridManager.c_grid()); const Wells* wells = wellsManager.c_wells(); // There is 3 wells in total in the deck at the 3rd schedule step. @@ -220,8 +220,8 @@ BOOST_AUTO_TEST_CASE(WellsEqual) { Opm::EclipseState eclipseState(deck, parseContext); Opm::GridManager gridManager(eclipseState.getInputGrid()); - Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid(), NULL); - Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid()); + Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid()); BOOST_CHECK(wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells(),false)); BOOST_CHECK(!wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells(),false)); @@ -235,8 +235,8 @@ BOOST_AUTO_TEST_CASE(ControlsEqual) { Opm::EclipseState eclipseState(deck, parseContext); Opm::GridManager gridManager(eclipseState.getInputGrid()); - Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid(), NULL); - Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid()); + Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid()); BOOST_CHECK(well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false)); BOOST_CHECK(well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false)); @@ -257,7 +257,7 @@ BOOST_AUTO_TEST_CASE(WellShutOK) { Opm::EclipseState eclipseState(deck, parseContext); Opm::GridManager gridManager(eclipseState.getInputGrid()); - Opm::WellsManager wellsManager2(eclipseState, 2, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager2(eclipseState, 2, *gridManager.c_grid()); // Shut wells are not added to the deck. i.e number of wells should be 2-1 BOOST_CHECK(wellsManager2.c_wells()->number_of_wells == 1); @@ -273,7 +273,7 @@ BOOST_AUTO_TEST_CASE(WellSTOPOK) { Opm::EclipseState eclipseState(deck, parseContext); Opm::GridManager gridManager(eclipseState.getInputGrid()); - Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid(), NULL); + Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid()); const Wells* wells = wellsManager.c_wells(); const struct WellControls* ctrls0 = wells->ctrls[0]; diff --git a/tests/wells_manager_data.data b/tests/wells_manager_data.data index 1c10c4ef..3d07c0c0 100755 --- a/tests/wells_manager_data.data +++ b/tests/wells_manager_data.data @@ -21,6 +21,8 @@ DZV TOPS 100*10 / +PERMX +500*1.0 / SCHEDULE diff --git a/tests/wells_manager_data_wellSTOP.data b/tests/wells_manager_data_wellSTOP.data index 3b4a6432..d10970c1 100755 --- a/tests/wells_manager_data_wellSTOP.data +++ b/tests/wells_manager_data_wellSTOP.data @@ -22,6 +22,9 @@ DZV TOPS 100*10 / +PERMX +500*1.0 / + SCHEDULE WELSPECS diff --git a/tests/wells_stopped.data b/tests/wells_stopped.data index 3a39481b..cfcda249 100644 --- a/tests/wells_stopped.data +++ b/tests/wells_stopped.data @@ -26,6 +26,8 @@ TOPS 100*10 / +PERMX +100*1.0 / SCHEDULE