My Project
BlackoilAquiferModel_impl.hpp
1 /*
2  Copyright 2017 TNO - Heat Transfer & Fluid Dynamics, Modelling & Optimization of the Subsurface
3  Copyright 2017 Statoil ASA.
4 
5  This file is part of the Open Porous Media project (OPM).
6 
7  OPM is free software: you can redistribute it and/or modify
8  it under the terms of the GNU General Public License as published by
9  the Free Software Foundation, either version 3 of the License, or
10  (at your option) any later version.
11 
12  OPM is distributed in the hope that it will be useful,
13  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15  GNU General Public License for more details.
16 
17  You should have received a copy of the GNU General Public License
18  along with OPM. If not, see <http://www.gnu.org/licenses/>.
19 */
20 
21 #include <opm/grid/utility/cartesianToCompressed.hpp>
22 #include "BlackoilAquiferModel.hpp"
23 
24 namespace Opm
25 {
26 
27 template <typename TypeTag>
28 BlackoilAquiferModel<TypeTag>::BlackoilAquiferModel(Simulator& simulator)
29  : simulator_(simulator)
30 {
31  // Grid needs to support Facetag
32  using Grid = std::remove_const_t<std::remove_reference_t<decltype(simulator.vanguard().grid())>>;
33  static_assert(SupportsFaceTag<Grid>::value, "Grid has to support assumptions about face tag.");
34 
35  init();
36 }
37 
38 template <typename TypeTag>
39 void
40 BlackoilAquiferModel<TypeTag>::initialSolutionApplied()
41 {
42  if (aquiferCarterTracyActive()) {
43  for (auto& aquifer : aquifers_CarterTracy) {
44  aquifer.initialSolutionApplied();
45  }
46  }
47  if (aquiferFetkovichActive()) {
48  for (auto& aquifer : aquifers_Fetkovich) {
49  aquifer.initialSolutionApplied();
50  }
51  }
52 
53  if (this->aquiferNumericalActive()) {
54  for (auto& aquifer : this->aquifers_numerical) {
55  aquifer.initialSolutionApplied();
56  }
57  }
58 }
59 
60 template <typename TypeTag>
61 void
62 BlackoilAquiferModel<TypeTag>::initFromRestart(const data::Aquifers& aquiferSoln)
63 {
64  if (this->aquiferCarterTracyActive()) {
65  for (auto& aquifer : this->aquifers_CarterTracy) {
66  aquifer.initFromRestart(aquiferSoln);
67  }
68  }
69 
70  if (this->aquiferFetkovichActive()) {
71  for (auto& aquifer : this->aquifers_Fetkovich) {
72  aquifer.initFromRestart(aquiferSoln);
73  }
74  }
75 
76  if (this->aquiferNumericalActive()) {
77  for (auto& aquifer : this->aquifers_numerical) {
78  aquifer.initFromRestart(aquiferSoln);
79  }
80  }
81 }
82 
83 template <typename TypeTag>
84 void
85 BlackoilAquiferModel<TypeTag>::beginEpisode()
86 {}
87 
88 template <typename TypeTag>
89 void
90 BlackoilAquiferModel<TypeTag>::beginIteration()
91 {}
92 
93 template <typename TypeTag>
94 void
95 BlackoilAquiferModel<TypeTag>::beginTimeStep()
96 {
97  if (aquiferCarterTracyActive()) {
98  for (auto& aquifer : aquifers_CarterTracy) {
99  aquifer.beginTimeStep();
100  }
101  }
102  if (aquiferFetkovichActive()) {
103  for (auto& aquifer : aquifers_Fetkovich) {
104  aquifer.beginTimeStep();
105  }
106  }
107 }
108 
109 template <typename TypeTag>
110 template <class Context>
111 void
112 BlackoilAquiferModel<TypeTag>::addToSource(RateVector& rates,
113  const Context& context,
114  unsigned spaceIdx,
115  unsigned timeIdx) const
116 {
117  if (aquiferCarterTracyActive()) {
118  for (auto& aquifer : aquifers_CarterTracy) {
119  aquifer.addToSource(rates, context, spaceIdx, timeIdx);
120  }
121  }
122  if (aquiferFetkovichActive()) {
123  for (auto& aquifer : aquifers_Fetkovich) {
124  aquifer.addToSource(rates, context, spaceIdx, timeIdx);
125  }
126  }
127 }
128 
129 template <typename TypeTag>
130 void
131 BlackoilAquiferModel<TypeTag>::endIteration()
132 {}
133 
134 template <typename TypeTag>
135 void
136 BlackoilAquiferModel<TypeTag>::endTimeStep()
137 {
138  if (aquiferCarterTracyActive()) {
139  for (auto& aquifer : aquifers_CarterTracy) {
140  aquifer.endTimeStep();
141  }
142  }
143  if (aquiferFetkovichActive()) {
144  for (auto& aquifer : aquifers_Fetkovich) {
145  aquifer.endTimeStep();
146  }
147  }
148  if (aquiferNumericalActive()) {
149  for (auto& aquifer : this->aquifers_numerical) {
150  aquifer.endTimeStep();
151  this->simulator_.vanguard().grid().comm().barrier();
152  }
153  }
154 }
155 
156 template <typename TypeTag>
157 void
158 BlackoilAquiferModel<TypeTag>::endEpisode()
159 {}
160 
161 template <typename TypeTag>
162 template <class Restarter>
163 void
164 BlackoilAquiferModel<TypeTag>::serialize(Restarter& /* res */)
165 {
166  // TODO (?)
167  throw std::logic_error("BlackoilAquiferModel::serialize() is not yet implemented");
168 }
169 
170 template <typename TypeTag>
171 template <class Restarter>
172 void
173 BlackoilAquiferModel<TypeTag>::deserialize(Restarter& /* res */)
174 {
175  // TODO (?)
176  throw std::logic_error("BlackoilAquiferModel::deserialize() is not yet implemented");
177 }
178 
179 // Initialize the aquifers in the deck
180 template <typename TypeTag>
181 void
182 BlackoilAquiferModel<TypeTag>::init()
183 {
184  const auto& aquifer = this->simulator_.vanguard().eclState().aquifer();
185 
186  if (!aquifer.active()) {
187  return;
188  }
189 
190 
191  const auto& connections = aquifer.connections();
192  for (const auto& aq : aquifer.ct()) {
193  aquifers_CarterTracy.emplace_back(connections[aq.aquiferID],
194  this->simulator_, aq);
195  }
196 
197  for (const auto& aq : aquifer.fetp()) {
198  aquifers_Fetkovich.emplace_back(connections[aq.aquiferID],
199  this->simulator_, aq);
200  }
201 
202  if (aquifer.hasNumericalAquifer()) {
203  const auto& num_aquifers = aquifer.numericalAquifers().aquifers();
204  const auto& ugrid = simulator_.vanguard().grid();
205  const int number_of_cells = simulator_.gridView().size(0);
206  const int* global_cell = UgGridHelpers::globalCell(ugrid);
207  const std::unordered_map<int, int> cartesian_to_compressed = cartesianToCompressed(number_of_cells,
208  global_cell);
209  for ([[maybe_unused]]const auto& [id, aqu] : num_aquifers) {
210  this->aquifers_numerical.emplace_back(aqu,
211  cartesian_to_compressed, this->simulator_, global_cell);
212  }
213  }
214 }
215 
216 template <typename TypeTag>
217 bool
218 BlackoilAquiferModel<TypeTag>::aquiferCarterTracyActive() const
219 {
220  return !aquifers_CarterTracy.empty();
221 }
222 
223 template <typename TypeTag>
224 bool
225 BlackoilAquiferModel<TypeTag>::aquiferFetkovichActive() const
226 {
227  return !aquifers_Fetkovich.empty();
228 }
229 
230 template<typename TypeTag>
231 bool
232 BlackoilAquiferModel<TypeTag>::aquiferNumericalActive() const
233 {
234  return !(this->aquifers_numerical.empty());
235 }
236 
237 template<typename TypeTag>
238 data::Aquifers BlackoilAquiferModel<TypeTag>::aquiferData() const
239 {
240  data::Aquifers data;
241  if (this->aquiferCarterTracyActive()) {
242  for (const auto& aqu : this->aquifers_CarterTracy) {
243  data.insert_or_assign(aqu.aquiferID(), aqu.aquiferData());
244  }
245  }
246 
247  if (this->aquiferFetkovichActive()) {
248  for (const auto& aqu : this->aquifers_Fetkovich) {
249  data.insert_or_assign(aqu.aquiferID(), aqu.aquiferData());
250  }
251  }
252 
253  if (this->aquiferNumericalActive()) {
254  for (const auto& aqu : this->aquifers_numerical) {
255  data.insert_or_assign(aqu.aquiferID(), aqu.aquiferData());
256  }
257  }
258 
259  return data;
260 }
261 } // namespace Opm
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27