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