41 using RateVector = GetPropType<TypeTag, Properties::RateVector>;
42 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
43 using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
44 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
45 using BlackoilIndices = GetPropType<TypeTag, Properties::Indices>;
47 static constexpr int numEq = BlackoilIndices::numEq;
48 using Eval = DenseAd::Evaluation<double, numEq>;
51 const Simulator& simulator,
52 const SingleAquiferFlux& aquifer)
54 , connections_ (connections)
55 , aquifer_data_ (aquifer)
56 , connection_flux_ (connections_.size(), Eval{0})
58 this->total_face_area_ = this->initializeConnections();
64 result.cumulative_flux_ = 1.0;
71 void computeFaceAreaFraction(
const std::vector<double>& total_face_area)
override
73 assert (total_face_area.size() >=
static_cast<std::vector<double>::size_type
>(this->aquiferID()));
75 this->area_fraction_ = this->totalFaceArea()
76 / total_face_area[this->aquiferID() - 1];
79 double totalFaceArea()
const override
81 return this->total_face_area_;
84 void updateAquifer(
const SingleAquiferFlux& aquifer)
86 aquifer_data_ = aquifer;
89 void initFromRestart(
const data::Aquifers& aquiferSoln)
override
91 auto xaqPos = aquiferSoln.find(this->aquiferID());
92 if (xaqPos == aquiferSoln.end()) {
96 this->cumulative_flux_ = this->area_fraction_ * xaqPos->second.volume;
99 void initialSolutionApplied()
override
102 void beginTimeStep()
override
105 void endTimeStep()
override
107 this->flux_rate_ = this->totalFluxRate();
108 this->cumulative_flux_ +=
109 this->flux_rate_ * this->simulator_.timeStepSize();
112 data::AquiferData aquiferData()
const override
114 data::AquiferData data;
116 data.aquiferID = this->aquifer_data_.id;
120 data.fluxRate = this->totalFluxRate();
122 data.volume = this->cumulative_flux_;
125 data.initPressure = 0.0;
130 void addToSource(RateVector& rates,
131 const unsigned cellIdx,
132 [[maybe_unused]]
const unsigned timeIdx)
override
134 const int idx = this->cellToConnectionIdx_[cellIdx];
139 const auto& model = this->simulator_.model();
141 const auto fw = this->aquifer_data_.flux;
143 this->connection_flux_[idx] = fw * this->connections_[idx].effective_facearea;
145 rates[BlackoilIndices::conti0EqIdx + compIdx_()]
146 += this->connection_flux_[idx] / model.dofTotalVolume(cellIdx);
149 template<
class Serializer>
150 void serializeOp(Serializer& serializer)
152 serializer(cumulative_flux_);
157 return this->cumulative_flux_ == rhs.cumulative_flux_;
161 const std::vector<Aquancon::AquancCell>& connections_;
163 SingleAquiferFlux aquifer_data_;
164 std::vector<Eval> connection_flux_{};
165 std::vector<int> cellToConnectionIdx_{};
167 double cumulative_flux_{};
168 double total_face_area_{0.0};
169 double area_fraction_{1.0};
171 double initializeConnections()
173 auto connected_face_area = 0.0;
175 this->cellToConnectionIdx_
176 .resize(this->simulator_.gridView().size(0), -1);
178 for (std::size_t idx = 0; idx < this->connections_.size(); ++idx) {
179 const auto global_index = this->connections_[idx].global_index;
180 const int cell_index = this->simulator_.vanguard()
181 .compressedIndexForInterior(global_index);
183 if (cell_index < 0) {
187 this->cellToConnectionIdx_[cell_index] = idx;
189 connected_face_area += this->connections_[idx].effective_facearea;
196 return connected_face_area;
199 double computeFaceAreaFraction(
const double connected_face_area)
const
201 const auto tot_face_area = this->simulator_.vanguard()
202 .grid().comm().sum(connected_face_area);
204 return (tot_face_area > 0.0)
205 ? connected_face_area / tot_face_area
212 if (this->co2store_or_h2store_())
213 return FluidSystem::oilCompIdx;
215 return FluidSystem::waterCompIdx;
218 double totalFluxRate()
const
220 return std::accumulate(this->connection_flux_.begin(),
221 this->connection_flux_.end(), 0.0,
222 [](
const double rate,
const auto& q)
224 return rate + getValue(q);