56 enum { dimWorld = GridView::dimensionworld };
57 enum { numPhases = FluidSystem::numPhases };
58 static constexpr int numEq = BlackoilIndices::numEq;
60 using Eval = DenseAd::Evaluation<double, numEq>;
70 , cumulative_flux_(0.0)
71 , init_pressure_ (
aquifer.numCells(), 0.0)
73 this->cell_to_aquifer_cell_idx_.resize(this->ebos_simulator_.gridView().size(0), -1);
88 this->checkConnectsToReservoir();
96 result.cumulative_flux_ = 2.0;
97 result.init_pressure_ = {3.0, 4.0};
103 void initFromRestart(
const data::Aquifers&
aquiferSoln)
override
109 if (this->connects_to_reservoir_) {
110 this->cumulative_flux_ =
xaqPos->second.volume;
116 this->init_pressure_ =
aqData->initPressure;
119 this->solution_set_from_restart_ =
true;
122 void beginTimeStep()
override {}
123 void addToSource(RateVector&,
const unsigned,
const unsigned)
override {}
125 void endTimeStep()
override
127 this->pressure_ = this->calculateAquiferPressure();
128 this->flux_rate_ = this->calculateAquiferFluxRate();
129 this->cumulative_flux_ += this->flux_rate_ * this->ebos_simulator_.timeStepSize();
132 data::AquiferData aquiferData()
const override
134 data::AquiferData data;
135 data.aquiferID = this->aquiferID();
136 data.pressure = this->pressure_;
137 data.fluxRate = this->flux_rate_;
138 data.volume = this->cumulative_flux_;
141 aquNum->initPressure = this->init_pressure_;
146 void initialSolutionApplied()
override
148 if (this->solution_set_from_restart_) {
152 this->pressure_ = this->calculateAquiferPressure(this->init_pressure_);
153 this->flux_rate_ = 0.;
154 this->cumulative_flux_ = 0.;
157 void computeFaceAreaFraction(
const std::vector<double>& )
override
160 double totalFaceArea()
const override
165 template<
class Serializer>
176 return this->flux_rate_ == rhs.flux_rate_ &&
177 this->cumulative_flux_ == rhs.cumulative_flux_ &&
178 this->init_pressure_ == rhs.init_pressure_ &&
179 this->pressure_ == rhs.pressure_;
182 double cumulativeFlux()
const
184 return this->cumulative_flux_;
188 void checkConnectsToReservoir()
195 elem_ctx.updateStencil(elem);
197 const auto cell_index = elem_ctx
198 .globalSpaceIndex(0, 0);
200 return this->cell_to_aquifer_cell_idx_[cell_index] == 0;
204 &&
"Internal error locating numerical aquifer's connecting cell");
206 this->connects_to_reservoir_ =
207 elemIt->partitionType() == Dune::InteriorEntity;
210 double calculateAquiferPressure()
const
212 auto capture = std::vector<double>(this->init_pressure_.size(), 0.0);
213 return this->calculateAquiferPressure(
capture);
216 double calculateAquiferPressure(std::vector<double>&
cell_pressure)
const
222 const auto&
gridView = this->ebos_simulator_.gridView();
223 OPM_BEGIN_PARALLEL_TRY_CATCH();
228 const std::size_t cell_index =
elem_ctx.globalSpaceIndex(0, 0);
229 const int idx = this->cell_to_aquifer_cell_idx_[cell_index];
234 elem_ctx.updatePrimaryIntensiveQuantities(0);
235 const auto&
iq0 =
elem_ctx.intensiveQuantities(0, 0);
236 const auto&
fs =
iq0.fluidState();
252 OPM_END_PARALLEL_TRY_CATCH(
"AquiferNumerical::calculateAquiferPressure() failed: ", this->ebos_simulator_.vanguard().grid().comm());
253 const auto& comm = this->ebos_simulator_.vanguard().grid().comm();
263 template <
class ElemCtx>
271 double calculateAquiferFluxRate()
const
275 if (! this->connects_to_reservoir_) {
280 const auto&
gridView = this->ebos_simulator_.gridView();
285 const std::size_t cell_index =
elem_ctx.globalSpaceIndex(0, 0);
286 const int idx = this->cell_to_aquifer_cell_idx_[cell_index];
300 const std::size_t
i =
face.interiorIndex();
301 const std::size_t
j =
face.exteriorIndex();
304 const std::size_t
J =
stencil.globalSpaceIndex(
j);
310 if (this->cell_to_aquifer_cell_idx_[
J] > 0) {
313 elem_ctx.updateAllIntensiveQuantities();
314 elem_ctx.updateAllExtensiveQuantities();
332 double cumulative_flux_;
333 std::vector<double> init_pressure_{};
335 bool solution_set_from_restart_ {
false};
336 bool connects_to_reservoir_ {
false};
339 std::vector<int> cell_to_aquifer_cell_idx_;