22#ifndef OPM_AQUIFERANALYTICAL_HEADER_INCLUDED
23#define OPM_AQUIFERANALYTICAL_HEADER_INCLUDED
25#include <dune/grid/common/partitionset.hh>
27#include <opm/common/ErrorMacros.hpp>
29#include <opm/input/eclipse/EclipseState/Aquifer/Aquancon.hpp>
31#include <opm/material/common/MathToolbox.hpp>
32#include <opm/material/densead/Evaluation.hpp>
33#include <opm/material/densead/Math.hpp>
34#include <opm/material/fluidstates/BlackOilFluidState.hpp>
39#include <opm/output/data/Aquifer.hpp>
54template <
typename TypeTag>
68 enum {
enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
69 enum {
enableBrine = getPropValue<TypeTag, Properties::EnableBrine>() };
70 enum {
enableVapwat = getPropValue<TypeTag, Properties::EnableVapwat>() };
75 static constexpr int numEq = BlackoilIndices::numEq;
83 BlackoilIndices::gasEnabled,
88 BlackoilIndices::numPhases>;
92 const std::vector<Aquancon::AquancCell>& connections,
106 assert (total_face_area.size() >=
static_cast<typename std::vector<Scalar>::size_type
>(this->aquiferID()));
108 const auto tfa = total_face_area[this->
aquiferID() - 1];
109 const auto eps_sqrt = std::sqrt(std::numeric_limits<Scalar>::epsilon());
111 if (tfa < eps_sqrt) {
116 this->faceArea_connected_.end(),
117 this->alphai_.begin(),
134 auto xaqPos = aquiferSoln.find(this->
aquiferID());
135 if (xaqPos == aquiferSoln.end())
141 this->
pa0_ = xaqPos->second.initPressure;
156 for (
const auto& elem : elements(this->
simulator_.gridView())) {
157 elemCtx.updatePrimaryStencil(elem);
159 const int cellIdx = elemCtx.globalSpaceIndex(0, 0);
164 elemCtx.updateIntensiveQuantities(0);
165 const auto& iq = elemCtx.intensiveQuantities(0, 0);
174 const unsigned cellIdx,
175 const unsigned timeIdx)
override
183 const auto& intQuants = model.intensiveQuantities(cellIdx, timeIdx);
189 rates[BlackoilIndices::conti0EqIdx +
compIdx_()]
190 += this->
Qai_[idx] / model.dofTotalVolume(cellIdx);
193 auto fs = intQuants.fluidState();
194 if (this->
Ta0_.has_value() && this->
Qai_[idx] > 0)
196 fs.setTemperature(this->
Ta0_.value());
197 typedef typename std::decay<
decltype(fs)>::type::Scalar FsScalar;
198 typename FluidSystem::template ParameterCache<FsScalar> paramCache;
199 const unsigned pvtRegionIdx = intQuants.pvtRegionIndex();
200 paramCache.setRegionIndex(pvtRegionIdx);
201 paramCache.setMaxOilSat(this->
simulator_.problem().maxOilSaturation(cellIdx));
202 paramCache.updatePhase(fs, this->
phaseIdx_());
203 const auto& h = FluidSystem::enthalpy(fs, paramCache, this->
phaseIdx_());
206 rates[BlackoilIndices::contiEnergyEqIdx]
207 += this->
Qai_[idx] *fs.enthalpy(this->
phaseIdx_()) * FluidSystem::referenceDensity( this->
phaseIdx_(), intQuants.pvtRegionIndex()) / model.dofTotalVolume(cellIdx);
217 template<
class Serializer>
245 return this->
simulator_.problem().gravity()[2];
251 return FluidSystem::oilCompIdx;
253 return FluidSystem::waterCompIdx;
269 this->
Qai_.resize(this->
size(), Scalar{0});
276 const auto& fs = intQuants.fluidState();
277 pressure_water.at(idx) = fs.pressure(this->
phaseIdx_());
284 const auto& fs = intQuants.fluidState();
285 pressure_water.at(idx) = fs.pressure(this->
phaseIdx_()).value();
296 const auto& gridView = this->
simulator_.vanguard().gridView();
297 for (std::size_t idx = 0; idx < this->
size(); ++idx) {
298 const auto global_index = this->
connections_[idx].global_index;
299 const int cell_index = this->
simulator_.vanguard().compressedIndex(global_index);
300 if (cell_index < 0) {
304 auto elemIt = gridView.template begin< 0>();
305 std::advance(elemIt, cell_index);
308 if (elemIt->partitionType() != Dune::InteriorEntity) {
316 FaceDir::DirEnum faceDirection;
319 const auto& elemMapper = this->
simulator_.model().dofMapper();
320 for (
const auto& elem : elements(gridView)) {
321 const unsigned cell_index = elemMapper.index(elem);
329 for (
const auto& intersection : intersections(gridView, elem)) {
331 if (! intersection.boundary()) {
335 switch (intersection.indexInInside()) {
337 faceDirection = FaceDir::XMinus;
340 faceDirection = FaceDir::XPlus;
343 faceDirection = FaceDir::YMinus;
346 faceDirection = FaceDir::YPlus;
349 faceDirection = FaceDir::ZMinus;
352 faceDirection = FaceDir::ZPlus;
355 OPM_THROW(std::logic_error,
356 "Internal error in initialization of aquifer.");
359 if (faceDirection == this->
connections_[idx].face_dir) {
373 const auto& gridView = this->
simulator_.vanguard().gridView();
374 for (std::size_t idx = 0; idx < this->
size(); ++idx) {
375 const int cell_index = this->
simulator_.vanguard()
377 if (cell_index < 0) {
381 auto elemIt = gridView.template begin< 0>();
382 std::advance(elemIt, cell_index);
385 if (elemIt->partitionType() != Dune::InteriorEntity) {
397 std::vector<Scalar> pw_aquifer;
398 Scalar water_pressure_reservoir;
401 const auto& gridView = this->
simulator_.gridView();
402 for (
const auto& elem : elements(gridView)) {
403 elemCtx.updatePrimaryStencil(elem);
405 const auto cellIdx = elemCtx.globalSpaceIndex(0, 0);
410 elemCtx.updatePrimaryIntensiveQuantities(0);
411 const auto& iq0 = elemCtx.intensiveQuantities(0, 0);
412 const auto& fs = iq0.fluidState();
414 water_pressure_reservoir = fs.pressure(this->
phaseIdx_()).value();
415 const auto water_density = fs.density(this->
phaseIdx_());
420 pw_aquifer.push_back(this->
alphai_[idx] *
421 (water_pressure_reservoir - water_density.value()*gdz));
425 const auto& comm = this->
simulator_.vanguard().grid().comm();
428 vals[0] = std::accumulate(this->
alphai_.begin(), this->alphai_.end(),
Scalar{0});
429 vals[1] = std::accumulate(pw_aquifer.begin(), pw_aquifer.end(),
Scalar{0});
433 return vals[1] / vals[0];
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:192
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:158
Defines a type tags and some fundamental properties all models.
Declares the properties required by the black oil model.
Definition: AquiferAnalytical.hpp:56
void updateCellPressure(std::vector< Scalar > &pressure_water, const int idx, const IntensiveQuantities &intQuants)
Definition: AquiferAnalytical.hpp:280
Scalar total_face_area_
Definition: AquiferAnalytical.hpp:454
Scalar area_fraction_
Definition: AquiferAnalytical.hpp:455
Scalar gravity_() const
Definition: AquiferAnalytical.hpp:243
std::vector< int > cellToConnectionIdx_
Definition: AquiferAnalytical.hpp:440
Scalar Tc_
Definition: AquiferAnalytical.hpp:449
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: AquiferAnalytical.hpp:58
void addToSource(RateVector &rates, const unsigned cellIdx, const unsigned timeIdx) override
Definition: AquiferAnalytical.hpp:173
virtual void calculateAquiferCondition()=0
Scalar totalFaceArea() const override
Definition: AquiferAnalytical.hpp:127
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: AquiferAnalytical.hpp:64
void computeFaceAreaFraction(const std::vector< Scalar > &total_face_area) override
Definition: AquiferAnalytical.hpp:104
int compIdx_() const
Definition: AquiferAnalytical.hpp:248
bool solution_set_from_restart_
Definition: AquiferAnalytical.hpp:459
void initializeConnectionMappings()
Definition: AquiferAnalytical.hpp:288
virtual Scalar aquiferDepth() const =0
std::size_t size() const
Definition: AquiferAnalytical.hpp:212
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: AquiferAnalytical.hpp:60
GetPropType< TypeTag, Properties::RateVector > RateVector
Definition: AquiferAnalytical.hpp:63
void initFromRestart(const data::Aquifers &aquiferSoln) override
Definition: AquiferAnalytical.hpp:132
std::vector< Scalar > faceArea_connected_
Definition: AquiferAnalytical.hpp:439
void initialSolutionApplied() override
Definition: AquiferAnalytical.hpp:146
void beginTimeStep() override
Definition: AquiferAnalytical.hpp:151
@ enableEnergy
Definition: AquiferAnalytical.hpp:68
std::vector< Scalar > alphai_
Definition: AquiferAnalytical.hpp:447
Scalar pa0_
Definition: AquiferAnalytical.hpp:450
DenseAd::Evaluation< Scalar, numEq > Eval
Definition: AquiferAnalytical.hpp:77
GetPropType< TypeTag, Properties::ElementMapper > ElementMapper
Definition: AquiferAnalytical.hpp:65
virtual void calculateInflowRate(int idx, const Simulator &simulator)=0
std::vector< Scalar > pressure_previous_
Definition: AquiferAnalytical.hpp:444
AquiferAnalytical(const int aqID, const std::vector< Aquancon::AquancCell > &connections, const Simulator &simulator)
Definition: AquiferAnalytical.hpp:91
void initQuantities()
Definition: AquiferAnalytical.hpp:256
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: AquiferAnalytical.hpp:59
Eval W_flux_
Definition: AquiferAnalytical.hpp:457
@ enableVapwat
Definition: AquiferAnalytical.hpp:70
virtual void calculateAquiferConstants()=0
std::vector< Eval > pressure_current_
Definition: AquiferAnalytical.hpp:445
@ enableBrine
Definition: AquiferAnalytical.hpp:69
@ enableSaltPrecipitation
Definition: AquiferAnalytical.hpp:73
static constexpr int numEq
Definition: AquiferAnalytical.hpp:75
std::optional< Scalar > Ta0_
Definition: AquiferAnalytical.hpp:451
std::vector< Scalar > cell_depth_
Definition: AquiferAnalytical.hpp:443
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: AquiferAnalytical.hpp:61
void serializeOp(Serializer &serializer)
Definition: AquiferAnalytical.hpp:218
bool has_active_connection_on_proc_
Definition: AquiferAnalytical.hpp:460
@ enableTemperature
Definition: AquiferAnalytical.hpp:67
const std::vector< Aquancon::AquancCell > connections_
Definition: AquiferAnalytical.hpp:436
std::vector< Eval > Qai_
Definition: AquiferAnalytical.hpp:446
@ has_disgas_in_water
Definition: AquiferAnalytical.hpp:71
Scalar calculateReservoirEquilibrium()
Definition: AquiferAnalytical.hpp:394
virtual ~AquiferAnalytical()
Definition: AquiferAnalytical.hpp:101
virtual void assignRestartData(const data::AquiferData &xaq)=0
void initializeConnectionDepths()
Definition: AquiferAnalytical.hpp:369
BlackOilFluidState< Eval, FluidSystem, enableTemperature, enableEnergy, BlackoilIndices::gasEnabled, enableVapwat, enableBrine, enableSaltPrecipitation, has_disgas_in_water, BlackoilIndices::numPhases > FluidState
Definition: AquiferAnalytical.hpp:88
void updateCellPressure(std::vector< Eval > &pressure_water, const int idx, const IntensiveQuantities &intQuants)
Definition: AquiferAnalytical.hpp:272
GetPropType< TypeTag, Properties::Indices > BlackoilIndices
Definition: AquiferAnalytical.hpp:62
bool operator==(const AquiferAnalytical &rhs) const
Definition: AquiferAnalytical.hpp:227
Scalar rhow_
Definition: AquiferAnalytical.hpp:452
Definition: AquiferInterface.hpp:35
const Simulator & simulator_
Definition: AquiferInterface.hpp:98
int phaseIdx_() const
Definition: AquiferInterface.hpp:88
bool co2store_or_h2store_() const
Definition: AquiferInterface.hpp:82
int aquiferID() const
Definition: AquiferInterface.hpp:79
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235