diff --git a/opm/simulators/aquifers/AquiferCarterTracy.hpp b/opm/simulators/aquifers/AquiferCarterTracy.hpp index 98168924a..8e96ca156 100644 --- a/opm/simulators/aquifers/AquiferCarterTracy.hpp +++ b/opm/simulators/aquifers/AquiferCarterTracy.hpp @@ -49,19 +49,19 @@ public: using Base::waterCompIdx; using Base::waterPhaseIdx; - AquiferCarterTracy(const Aquancon::AquanconOutput& connection, + AquiferCarterTracy(const std::vector& connections, const std::unordered_map& cartesian_to_compressed, const Simulator& ebosSimulator, const AquiferCT::AQUCT_data& aquct_data) - : Base(connection, cartesian_to_compressed, ebosSimulator) + : Base(aquct_data.aquiferID, connections, cartesian_to_compressed, ebosSimulator) , aquct_data_(aquct_data) { } void endTimeStep() override { - for (const auto& Qai : Base::Qai_) { - Base::W_flux_ += Qai * Base::ebos_simulator_.timeStepSize(); + for (const auto& q : this->Qai_) { + this->W_flux_ += q * this->ebos_simulator_.timeStepSize(); } } @@ -75,17 +75,16 @@ protected: // This function is used to initialize and calculate the alpha_i for each grid connection to the aquifer inline void initializeConnections() override { - const auto& eclState = Base::ebos_simulator_.vanguard().eclState(); - const auto& ugrid = Base::ebos_simulator_.vanguard().grid(); + const auto& eclState = this->ebos_simulator_.vanguard().eclState(); + const auto& ugrid = this->ebos_simulator_.vanguard().grid(); const auto& grid = eclState.getInputGrid(); - Base::cell_idx_ = this->connection_.global_index; auto globalCellIdx = ugrid.globalCell(); // We hack the cell depth values for now. We can actually get it from elementcontext pos - Base::cell_depth_.resize(Base::cell_idx_.size(), aquct_data_.d0); - Base::alphai_.resize(Base::cell_idx_.size(), 1.0); - Base::faceArea_connected_.resize(Base::cell_idx_.size(), 0.0); + this->cell_depth_.resize(this->size(), aquct_data_.d0); + this->alphai_.resize(this->size(), 1.0); + this->faceArea_connected_.resize(this->size(), 0.0); auto cell2Faces = Opm::UgGridHelpers::cell2Faces(ugrid); auto faceCells = Opm::UgGridHelpers::faceCells(ugrid); @@ -95,10 +94,10 @@ protected: // denom_face_areas is the sum of the areas connected to an aquifer Scalar denom_face_areas = 0.; - Base::cellToConnectionIdx_.resize(Base::ebos_simulator_.gridView().size(/*codim=*/0), -1); - for (size_t idx = 0; idx < Base::cell_idx_.size(); ++idx) { - const int cell_index = Base::cartesian_to_compressed_.at(Base::cell_idx_[idx]); - Base::cellToConnectionIdx_[cell_index] = idx; + this->cellToConnectionIdx_.resize(this->ebos_simulator_.gridView().size(/*codim=*/0), -1); + for (size_t idx = 0; idx < this->size(); ++idx) { + const int cell_index = this->cartesian_to_compressed_.at(this->connections_[idx].global_index); + this->cellToConnectionIdx_[cell_index] = idx; const auto cellFacesRange = cell2Faces[cell_index]; for (auto cellFaceIter = cellFacesRange.begin(); cellFaceIter != cellFacesRange.end(); ++cellFaceIter) { @@ -132,21 +131,21 @@ protected: "Initialization of Aquifer Carter Tracy problem. Make sure faceTag is correctly defined"); } - if (faceDirection == this->connection_.reservoir_face_dir.at(idx)) { - Base::faceArea_connected_.at(idx) = Base::getFaceArea(faceCells, ugrid, faceIdx, idx); - denom_face_areas += (this->connection_.influx_multiplier.at(idx) * Base::faceArea_connected_.at(idx)); + if (faceDirection == this->connections_[idx].face_dir) { + this->faceArea_connected_.at(idx) = this->getFaceArea(faceCells, ugrid, faceIdx, idx); + denom_face_areas += (this->connections_[idx].influx_mult * this->faceArea_connected_.at(idx)); } } - auto cellCenter = grid.getCellCenter(Base::cell_idx_.at(idx)); - Base::cell_depth_.at(idx) = cellCenter[2]; + auto cellCenter = grid.getCellCenter(this->connections_[idx].global_index); + this->cell_depth_.at(idx) = cellCenter[2]; } const double eps_sqrt = std::sqrt(std::numeric_limits::epsilon()); - for (size_t idx = 0; idx < Base::cell_idx_.size(); ++idx) { - Base::alphai_.at(idx) = (denom_face_areas < eps_sqrt) + for (size_t idx = 0; idx < this->size(); ++idx) { + this->alphai_.at(idx) = (denom_face_areas < eps_sqrt) ? // Prevent no connection NaNs due to division by zero 0. - : (this->connection_.influx_multiplier.at(idx) * Base::faceArea_connected_.at(idx)) / denom_face_areas; + : (this->connections_[idx].influx_mult * this->faceArea_connected_.at(idx)) / denom_face_areas; } } @@ -165,22 +164,22 @@ protected: inline Scalar dpai(int idx) { - Scalar dp = Base::pa0_ - + Base::rhow_.at(idx).value() * Base::gravity_() * (Base::cell_depth_.at(idx) - aquct_data_.d0) - - Base::pressure_previous_.at(idx); + Scalar dp = this->pa0_ + + this->rhow_.at(idx).value() * this->gravity_() * (this->cell_depth_.at(idx) - aquct_data_.d0) + - this->pressure_previous_.at(idx); return dp; } // This function implements Eqs 5.8 and 5.9 of the EclipseTechnicalDescription inline void calculateEqnConstants(Scalar& a, Scalar& b, const int idx, const Simulator& simulator) { - const Scalar td_plus_dt = (simulator.timeStepSize() + simulator.time()) / Base::Tc_; - const Scalar td = simulator.time() / Base::Tc_; + const Scalar td_plus_dt = (simulator.timeStepSize() + simulator.time()) / this->Tc_; + const Scalar td = simulator.time() / this->Tc_; Scalar PItdprime = 0.; Scalar PItd = 0.; getInfluenceTableValues(PItd, PItdprime, td_plus_dt); - a = 1.0 / Base::Tc_ * ((beta_ * dpai(idx)) - (Base::W_flux_.value() * PItdprime)) / (PItd - td * PItdprime); - b = beta_ / (Base::Tc_ * (PItd - td * PItdprime)); + a = 1.0 / this->Tc_ * ((beta_ * dpai(idx)) - (this->W_flux_.value() * PItdprime)) / (PItd - td * PItdprime); + b = beta_ / (this->Tc_ * (PItd - td * PItdprime)); } // This function implements Eq 5.7 of the EclipseTechnicalDescription @@ -188,8 +187,8 @@ protected: { Scalar a, b; calculateEqnConstants(a, b, idx, simulator); - Base::Qai_.at(idx) - = Base::alphai_.at(idx) * (a - b * (Base::pressure_current_.at(idx) - Base::pressure_previous_.at(idx))); + this->Qai_.at(idx) + = this->alphai_.at(idx) * (a - b * (this->pressure_current_.at(idx) - this->pressure_previous_.at(idx))); } inline void calculateAquiferConstants() override @@ -198,7 +197,7 @@ protected: beta_ = aquct_data_.c2 * aquct_data_.h * aquct_data_.theta * aquct_data_.phi_aq * aquct_data_.C_t * aquct_data_.r_o * aquct_data_.r_o; // We calculate the time constant - Base::Tc_ = mu_w_ * aquct_data_.phi_aq * aquct_data_.C_t * aquct_data_.r_o * aquct_data_.r_o + this->Tc_ = mu_w_ * aquct_data_.phi_aq * aquct_data_.C_t * aquct_data_.r_o * aquct_data_.r_o / (aquct_data_.k_a * aquct_data_.c1); } @@ -206,17 +205,17 @@ protected: { int pvttableIdx = aquct_data_.pvttableID - 1; - Base::rhow_.resize(Base::cell_idx_.size(), 0.); + this->rhow_.resize(this->size(), 0.); if (!aquct_data_.p0.first) { - Base::pa0_ = calculateReservoirEquilibrium(); + this->pa0_ = calculateReservoirEquilibrium(); } else { - Base::pa0_ = aquct_data_.p0.second; + this->pa0_ = aquct_data_.p0.second; } // use the thermodynamic state of the first active cell as a // reference. there might be better ways to do this... - ElementContext elemCtx(Base::ebos_simulator_); - auto elemIt = Base::ebos_simulator_.gridView().template begin(); + ElementContext elemCtx(this->ebos_simulator_); + auto elemIt = this->ebos_simulator_.gridView().template begin(); elemCtx.updatePrimaryStencil(*elemIt); elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0); const auto& iq0 = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0); @@ -227,7 +226,7 @@ protected: fs_aquifer.assign(iq0.fluidState()); Eval temperature_aq, pa0_mean; temperature_aq = fs_aquifer.temperature(0); - pa0_mean = Base::pa0_; + pa0_mean = this->pa0_; Eval mu_w_aquifer = FluidSystem::waterPvt().viscosity(pvttableIdx, temperature_aq, pa0_mean); mu_w_ = mu_w_aquifer.value(); } @@ -240,8 +239,8 @@ protected: std::vector pw_aquifer; Scalar water_pressure_reservoir; - ElementContext elemCtx(Base::ebos_simulator_); - const auto& gridView = Base::ebos_simulator_.gridView(); + ElementContext elemCtx(this->ebos_simulator_); + const auto& gridView = this->ebos_simulator_.gridView(); auto elemIt = gridView.template begin(); const auto& elemEndIt = gridView.template end(); for (; elemIt != elemEndIt; ++elemIt) { @@ -249,7 +248,7 @@ protected: elemCtx.updatePrimaryStencil(elem); size_t cellIdx = elemCtx.globalSpaceIndex(/*spaceIdx=*/0, /*timeIdx=*/0); - int idx = Base::cellToConnectionIdx_[cellIdx]; + int idx = this->cellToConnectionIdx_[cellIdx]; if (idx < 0) continue; @@ -258,11 +257,11 @@ protected: const auto& fs = iq0.fluidState(); water_pressure_reservoir = fs.pressure(waterPhaseIdx).value(); - Base::rhow_[idx] = fs.density(waterPhaseIdx); + this->rhow_[idx] = fs.density(waterPhaseIdx); pw_aquifer.push_back( (water_pressure_reservoir - - Base::rhow_[idx].value() * Base::gravity_() * (Base::cell_depth_[idx] - aquct_data_.d0)) - * Base::alphai_[idx]); + - this->rhow_[idx].value() * this->gravity_() * (this->cell_depth_[idx] - aquct_data_.d0)) + * this->alphai_[idx]); } // We take the average of the calculated equilibrium pressures. diff --git a/opm/simulators/aquifers/AquiferFetkovich.hpp b/opm/simulators/aquifers/AquiferFetkovich.hpp index aacb290d4..b34d437c7 100644 --- a/opm/simulators/aquifers/AquiferFetkovich.hpp +++ b/opm/simulators/aquifers/AquiferFetkovich.hpp @@ -51,19 +51,19 @@ public: using Base::waterCompIdx; using Base::waterPhaseIdx; - AquiferFetkovich(const Aquancon::AquanconOutput& connection, + AquiferFetkovich(const std::vector& connections, const std::unordered_map& cartesian_to_compressed, const Simulator& ebosSimulator, const Aquifetp::AQUFETP_data& aqufetp_data) - : Base(connection, cartesian_to_compressed, ebosSimulator) + : Base(aqufetp_data.aquiferID, connections, cartesian_to_compressed, ebosSimulator) , aqufetp_data_(aqufetp_data) { } void endTimeStep() override { - for (const auto& Qai : Base::Qai_) { - Base::W_flux_ += Qai * Base::ebos_simulator_.timeStepSize(); + for (const auto& q : this->Qai_) { + this->W_flux_ += q * this->ebos_simulator_.timeStepSize(); aquifer_pressure_ = aquiferPressure(); } } @@ -76,18 +76,17 @@ protected: inline void initializeConnections() override { - const auto& eclState = Base::ebos_simulator_.vanguard().eclState(); - const auto& ugrid = Base::ebos_simulator_.vanguard().grid(); + const auto& eclState = this->ebos_simulator_.vanguard().eclState(); + const auto& ugrid = this->ebos_simulator_.vanguard().grid(); const auto& grid = eclState.getInputGrid(); - Base::cell_idx_ = this->connection_.global_index; auto globalCellIdx = ugrid.globalCell(); // We hack the cell depth values for now. We can actually get it from elementcontext pos - Base::cell_depth_.resize(Base::cell_idx_.size(), aqufetp_data_.d0); - Base::alphai_.resize(Base::cell_idx_.size(), 1.0); - Base::faceArea_connected_.resize(Base::cell_idx_.size(), 0.0); + this->cell_depth_.resize(this->size(), aqufetp_data_.d0); + this->alphai_.resize(this->size(), 1.0); + this->faceArea_connected_.resize(this->size(), 0.0); auto cell2Faces = Opm::UgGridHelpers::cell2Faces(ugrid); auto faceCells = Opm::UgGridHelpers::faceCells(ugrid); @@ -97,14 +96,16 @@ protected: // denom_face_areas is the sum of the areas connected to an aquifer Scalar denom_face_areas = 0.; - Base::cellToConnectionIdx_.resize(Base::ebos_simulator_.gridView().size(/*codim=*/0), -1); - for (size_t idx = 0; idx < Base::cell_idx_.size(); ++idx) { - const int cell_index = Base::cartesian_to_compressed_.at(Base::cell_idx_[idx]); - Base::cellToConnectionIdx_[cell_index] = idx; - const auto cellCenter = grid.getCellCenter(Base::cell_idx_.at(idx)); - Base::cell_depth_.at(idx) = cellCenter[2]; + this->cellToConnectionIdx_.resize(this->ebos_simulator_.gridView().size(/*codim=*/0), -1); + for (size_t idx = 0; idx < this->size(); ++idx) { + const auto global_index = this->connections_[idx].global_index; + const int cell_index = this->cartesian_to_compressed_.at(global_index); - if (!this->connection_.influx_coeff[idx]) { // influx_coeff is defaulted + this->cellToConnectionIdx_[cell_index] = idx; + const auto cellCenter = grid.getCellCenter(global_index); + this->cell_depth_.at(idx) = cellCenter[2]; + + if (!this->connections_[idx].influx_coeff.first) { // influx_coeff is defaulted const auto cellFacesRange = cell2Faces[cell_index]; for (auto cellFaceIter = cellFacesRange.begin(); cellFaceIter != cellFacesRange.end(); ++cellFaceIter) { // The index of the face in the compressed grid @@ -137,24 +138,23 @@ protected: "Initialization of Aquifer problem. Make sure faceTag is correctly defined"); } - if (faceDirection == this->connection_.reservoir_face_dir.at(idx)) { - Base::faceArea_connected_.at(idx) - = Base::getFaceArea(faceCells, ugrid, faceIdx, idx); + if (faceDirection == this->connections_[idx].face_dir) { + this->faceArea_connected_[idx] = this->getFaceArea(faceCells, ugrid, faceIdx, idx); break; } } } else { - Base::faceArea_connected_.at(idx) = *this->connection_.influx_coeff[idx]; + this->faceArea_connected_.at(idx) = this->connections_[idx].influx_coeff.second; } - denom_face_areas += (this->connection_.influx_multiplier.at(idx) * Base::faceArea_connected_.at(idx)); + denom_face_areas += (this->connections_[idx].influx_mult * this->faceArea_connected_.at(idx)); } const double eps_sqrt = std::sqrt(std::numeric_limits::epsilon()); - for (size_t idx = 0; idx < Base::cell_idx_.size(); ++idx) { - Base::alphai_.at(idx) = (denom_face_areas < eps_sqrt) + for (size_t idx = 0; idx < this->size(); ++idx) { + this->alphai_.at(idx) = (denom_face_areas < eps_sqrt) ? // Prevent no connection NaNs due to division by zero 0. - : (this->connection_.influx_multiplier.at(idx) * Base::faceArea_connected_.at(idx)) / denom_face_areas; + : (this->connections_[idx].influx_mult * this->faceArea_connected_.at(idx)) / denom_face_areas; } } @@ -170,45 +170,45 @@ protected: inline Eval dpai(int idx) { - const Eval dp = aquifer_pressure_ - Base::pressure_current_.at(idx) - + Base::rhow_[idx] * Base::gravity_() * (Base::cell_depth_[idx] - aqufetp_data_.d0); + const Eval dp = aquifer_pressure_ - this->pressure_current_.at(idx) + + this->rhow_[idx] * this->gravity_() * (this->cell_depth_[idx] - aqufetp_data_.d0); return dp; } // This function implements Eq 5.12 of the EclipseTechnicalDescription inline Scalar aquiferPressure() { - Scalar Flux = Base::W_flux_.value(); - Scalar pa_ = Base::pa0_ - Flux / (aqufetp_data_.C_t * aqufetp_data_.V0); + Scalar Flux = this->W_flux_.value(); + Scalar pa_ = this->pa0_ - Flux / (aqufetp_data_.C_t * aqufetp_data_.V0); return pa_; } inline void calculateAquiferConstants() override { - Base::Tc_ = (aqufetp_data_.C_t * aqufetp_data_.V0) / aqufetp_data_.J; + this->Tc_ = (aqufetp_data_.C_t * aqufetp_data_.V0) / aqufetp_data_.J; } // This function implements Eq 5.14 of the EclipseTechnicalDescription inline void calculateInflowRate(int idx, const Simulator& simulator) override { - const Scalar td_Tc_ = simulator.timeStepSize() / Base::Tc_; + const Scalar td_Tc_ = simulator.timeStepSize() / this->Tc_; const Scalar coef = (1 - exp(-td_Tc_)) / td_Tc_; - Base::Qai_.at(idx) = Base::alphai_[idx] * aqufetp_data_.J * dpai(idx) * coef; + this->Qai_.at(idx) = this->alphai_[idx] * aqufetp_data_.J * dpai(idx) * coef; } inline void calculateAquiferCondition() override { - Base::rhow_.resize(Base::cell_idx_.size(), 0.); + this->rhow_.resize(this->size(), 0.); if (this->solution_set_from_restart_) { return; } if (!aqufetp_data_.p0.first) { - Base::pa0_ = calculateReservoirEquilibrium(); + this->pa0_ = calculateReservoirEquilibrium(); } else { - Base::pa0_ = aqufetp_data_.p0.second; + this->pa0_ = aqufetp_data_.p0.second; } - aquifer_pressure_ = Base::pa0_; + aquifer_pressure_ = this->pa0_; } inline Scalar calculateReservoirEquilibrium() override @@ -217,15 +217,15 @@ protected: std::vector pw_aquifer; Scalar water_pressure_reservoir; - ElementContext elemCtx(Base::ebos_simulator_); - const auto& gridView = Base::ebos_simulator_.gridView(); + ElementContext elemCtx(this->ebos_simulator_); + const auto& gridView = this->ebos_simulator_.gridView(); auto elemIt = gridView.template begin(); const auto& elemEndIt = gridView.template end(); for (; elemIt != elemEndIt; ++elemIt) { const auto& elem = *elemIt; elemCtx.updatePrimaryStencil(elem); size_t cellIdx = elemCtx.globalSpaceIndex(/*spaceIdx=*/0, /*timeIdx=*/0); - int idx = Base::cellToConnectionIdx_[cellIdx]; + int idx = this->cellToConnectionIdx_[cellIdx]; if (idx < 0) continue; @@ -234,11 +234,11 @@ protected: const auto& fs = iq0.fluidState(); water_pressure_reservoir = fs.pressure(waterPhaseIdx).value(); - Base::rhow_[idx] = fs.density(waterPhaseIdx); + this->rhow_[idx] = fs.density(waterPhaseIdx); pw_aquifer.push_back( (water_pressure_reservoir - - Base::rhow_[idx].value() * Base::gravity_() * (Base::cell_depth_[idx] - aqufetp_data_.d0)) - * Base::alphai_[idx]); + - this->rhow_[idx].value() * this->gravity_() * (this->cell_depth_[idx] - aqufetp_data_.d0)) + * this->alphai_[idx]); } // We take the average of the calculated equilibrium pressures. diff --git a/opm/simulators/aquifers/AquiferInterface.hpp b/opm/simulators/aquifers/AquiferInterface.hpp index 981364065..de3a5a3a7 100644 --- a/opm/simulators/aquifers/AquiferInterface.hpp +++ b/opm/simulators/aquifers/AquiferInterface.hpp @@ -73,16 +73,15 @@ public: static const auto waterPhaseIdx = FluidSystem::waterPhaseIdx; // Constructor - AquiferInterface(const Aquancon::AquanconOutput& connection, + AquiferInterface(int aqID, + const std::vector& connections, const std::unordered_map& cartesian_to_compressed, const Simulator& ebosSimulator) - : connection_(connection) + : aquiferID(aqID) + , connections_(connections) , ebos_simulator_(ebosSimulator) , cartesian_to_compressed_(cartesian_to_compressed) { - assert(this->connection_.influx_coeff.size() == this->connection_.global_index.size()); - assert(this->connection_.influx_coeff.size() == this->connection_.influx_multiplier.size()); - assert(this->connection_.influx_multiplier.size() == this->connection_.reservoir_face_dir.size()); } // Deconstructor @@ -94,16 +93,13 @@ public: { auto xaqPos = std::find_if(aquiferSoln.begin(), aquiferSoln.end(), [this](const data::AquiferData& xaq) -> bool { - return xaq.aquiferID == this->connection_.aquiferID; + return xaq.aquiferID == this->aquiferID; }); - if (xaqPos == aquiferSoln.end()) { - // No restart value applies to this aquifer. Nothing to do. + if (xaqPos == aquiferSoln.end()) return; - } this->assignRestartData(*xaqPos); - this->W_flux_ = xaqPos->volume; this->pa0_ = xaqPos->initPressure; this->solution_set_from_restart_ = true; @@ -155,6 +151,12 @@ public: += Qai_[idx] / context.dofVolume(spaceIdx, timeIdx); } + + std::size_t size() const { + return this->connections_.size(); + } + + protected: inline Scalar gravity_() const { @@ -174,9 +176,9 @@ protected: calculateAquiferCondition(); calculateAquiferConstants(); - pressure_previous_.resize(cell_idx_.size(), 0.); - pressure_current_.resize(cell_idx_.size(), 0.); - Qai_.resize(cell_idx_.size(), 0.0); + pressure_previous_.resize(this->connections_.size(), 0.); + pressure_current_.resize(this->connections_.size(), 0.); + Qai_.resize(this->connections_.size(), 0.0); } inline void @@ -214,7 +216,7 @@ protected: const auto cellNeighbour1 = faceCells(faceIdx, 1); const auto defaultFaceArea = Opm::UgGridHelpers::faceArea(ugrid, faceIdx); const auto calculatedFaceArea - = (!this->connection_.influx_coeff.at(idx)) ? defaultFaceArea : *(this->connection_.influx_coeff.at(idx)); + = (!this->connections_[idx].influx_coeff.first) ? defaultFaceArea : this->connections_[idx].influx_coeff.second; faceArea = (cellNeighbour0 * cellNeighbour1 > 0) ? 0. : calculatedFaceArea; if (cellNeighbour1 == 0) { faceArea = (cellNeighbour0 < 0) ? faceArea : 0.; @@ -226,12 +228,12 @@ protected: virtual void endTimeStep() = 0; - const Aquancon::AquanconOutput connection_; + const int aquiferID; + const std::vector connections_; const Simulator& ebos_simulator_; const std::unordered_map cartesian_to_compressed_; // Grid variables - std::vector cell_idx_; std::vector faceArea_connected_; std::vector cellToConnectionIdx_; // Quantities at each grid id diff --git a/opm/simulators/aquifers/BlackoilAquiferModel_impl.hpp b/opm/simulators/aquifers/BlackoilAquiferModel_impl.hpp index 58a82199f..aa677d74c 100644 --- a/opm/simulators/aquifers/BlackoilAquiferModel_impl.hpp +++ b/opm/simulators/aquifers/BlackoilAquiferModel_impl.hpp @@ -158,20 +158,15 @@ BlackoilAquiferModel::init() const AquiferCT aquiferct = AquiferCT(eclState.getTableManager(), deck); const Aquancon aquifer_connect = Aquancon(eclState.getInputGrid(), deck); - std::vector aquifersData = aquiferct.data(); - std::vector aquifer_connection = aquifer_connect.getAquOutput(); - - assert(aquifersData.size() == aquifer_connection.size()); const auto& ugrid = simulator_.vanguard().grid(); const auto& gridView = simulator_.gridView(); const int number_of_cells = gridView.size(0); cartesian_to_compressed_ = cartesianToCompressed(number_of_cells, Opm::UgGridHelpers::globalCell(ugrid)); - for (size_t i = 0; i < aquifersData.size(); ++i) { - aquifers_CarterTracy.push_back(AquiferCarterTracy( - aquifer_connection[i], cartesian_to_compressed_, this->simulator_, aquifersData[i])); - } + for (const auto& aquifer : aquiferct) + aquifers_CarterTracy.push_back(AquiferCarterTracy(aquifer_connect[aquifer.aquiferID], cartesian_to_compressed_, this->simulator_, aquifer)); + } if (comm.rank() == 0) has = deck.hasKeyword("AQUFETP"); @@ -188,20 +183,15 @@ BlackoilAquiferModel::init() const Aquifetp aquifetp = Aquifetp(deck); const Aquancon aquifer_connect = Aquancon(eclState.getInputGrid(), deck); - std::vector aquifersData = aquifetp.data(); - std::vector aquifer_connection = aquifer_connect.getAquOutput(); - - assert(aquifersData.size() == aquifer_connection.size()); const auto& ugrid = simulator_.vanguard().grid(); const auto& gridView = simulator_.gridView(); const int number_of_cells = gridView.size(0); cartesian_to_compressed_ = cartesianToCompressed(number_of_cells, Opm::UgGridHelpers::globalCell(ugrid)); - for (size_t i = 0; i < aquifersData.size(); ++i) { - aquifers_Fetkovich.push_back(AquiferFetkovich( - aquifer_connection[i], cartesian_to_compressed_, this->simulator_, aquifersData[i])); - } + for (const auto& aquifer : aquifetp) + aquifers_Fetkovich.push_back(AquiferFetkovich(aquifer_connect[aquifer.aquiferID], cartesian_to_compressed_, this->simulator_, aquifer)); + } } template