My Project
Loading...
Searching...
No Matches
FlowBaseVanguard.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_FLOW_BASE_VANGUARD_HPP
28#define OPM_FLOW_BASE_VANGUARD_HPP
29
30#include <dune/grid/common/partitionset.hh>
31
32#include <opm/grid/common/GridEnums.hpp>
33#include <opm/grid/common/CartesianIndexMapper.hpp>
34#include <opm/grid/LookUpCellCentroid.hh>
35
36#include <opm/input/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferCell.hpp>
37#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
38
44
45#include <opm/simulators/flow/BlackoilModelParameters.hpp>
47
48#include <array>
49#include <cstddef>
50#include <unordered_map>
51#include <vector>
52
53namespace Opm {
54template <class TypeTag>
55class FlowBaseVanguard;
56template<typename Grid, typename GridView> struct LookUpCellCentroid;
57}
58
59namespace Opm::Properties {
60
61namespace TTag {
62
64
65}
66
67// declare the properties required by the for the ecl simulator vanguard
68template<class TypeTag, class MyTypeTag>
69struct EquilGrid { using type = UndefinedProperty; };
70
71} // namespace Opm::Properties
72
73namespace Opm {
74
80template <class TypeTag>
81class FlowBaseVanguard : public BaseVanguard<TypeTag>,
83{
85 using Implementation = GetPropType<TypeTag, Properties::Vanguard>;
89
90 enum { enableExperiments = getPropValue<TypeTag, Properties::EnableExperiments>() };
91
92public:
95
96protected:
97 static const int dimension = Grid::dimension;
98 static const int dimensionworld = Grid::dimensionworld;
99 using Element = typename GridView::template Codim<0>::Entity;
101
102public:
106 static void registerParameters()
107 {
108 FlowGenericVanguard::registerParameters_<Scalar>();
109 }
110
117 explicit FlowBaseVanguard(Simulator& simulator)
118 : ParentType(simulator)
119 {
120#if HAVE_MPI
121 imbalanceTol_ = Parameters::Get<Parameters::ImbalanceTol<Scalar>>();
122
123 zoltanImbalanceTolSet_ = Parameters::IsSet<Parameters::ZoltanImbalanceTol<Scalar>>();
124 zoltanImbalanceTol_ = Parameters::Get<Parameters::ZoltanImbalanceTol<Scalar>>();
125#endif
126 enableExperiments_ = enableExperiments;
127
128 init();
129 }
130
131 const CartesianIndexMapper& cartesianMapper() const
132 { return asImp_().cartesianIndexMapper(); }
133
137 const std::array<int, dimension>& cartesianDimensions() const
138 { return asImp_().cartesianIndexMapper().cartesianDimensions(); }
139
143 int cartesianSize() const
144 { return asImp_().cartesianIndexMapper().cartesianSize(); }
145
150 { return asImp_().equilCartesianIndexMapper().cartesianSize(); }
151
155 unsigned cartesianIndex(unsigned compressedCellIdx) const
156 { return asImp_().cartesianIndexMapper().cartesianIndex(compressedCellIdx); }
157
161 unsigned cartesianIndex(const std::array<int,dimension>& coords) const
162 {
163 unsigned cartIndex = coords[0];
164 int factor = cartesianDimensions()[0];
165 for (unsigned i = 1; i < dimension; ++i) {
166 cartIndex += coords[i]*factor;
168 }
169
170 return cartIndex;
171 }
172
180 {
183 return index_pair->second;
184 else
185 return -1;
186 }
187
195 {
197 if (index_pair == cartesianToCompressed_.end() ||
198 !is_interior_[index_pair->second])
199 {
200 return -1;
201 }
202 else
203 {
204 return index_pair->second;
205 }
206 }
213 void cartesianCoordinate(unsigned cellIdx, std::array<int,3>& ijk) const
214 { return asImp_().cartesianIndexMapper().cartesianCoordinate(cellIdx, ijk); }
215
220 { return asImp_().equilCartesianIndexMapper().cartesianIndex(compressedEquilCellIdx); }
221
228 void equilCartesianCoordinate(unsigned cellIdx, std::array<int,3>& ijk) const
229 { return asImp_().equilCartesianIndexMapper().cartesianCoordinate(cellIdx, ijk); }
230
231
238 Scalar cellCenterDepth(unsigned globalSpaceIdx) const
239 {
241 }
242
243 const std::vector<Scalar>& cellCenterDepths() const
244 {
245 return cellCenterDepth_;
246 }
247
255 Scalar cellThickness(unsigned globalSpaceIdx) const
256 {
257 assert(!cellThickness_.empty());
259 }
260
266 std::size_t globalNumCells() const
267 {
268 const auto& grid = asImp_().grid();
269 if (grid.comm().size() == 1)
270 {
271 return grid.leafGridView().size(0);
272 }
273 const auto& gridView = grid.leafGridView();
274 constexpr int codim = 0;
275 constexpr auto Part = Dune::Interior_Partition;
276 auto local_cells = std::distance(gridView.template begin<codim, Part>(),
277 gridView.template end<codim, Part>());
278 return grid.comm().sum(local_cells);
279 }
280
281 void setupCartesianToCompressed_() {
282 this->updateCartesianToCompressedMapping_();
283 }
284
285protected:
295 template<class CartMapper>
296 std::function<std::array<double,dimensionworld>(int)>
297 cellCentroids_(const CartMapper& cartMapper, const bool& isCpGrid) const
298 {
299 return [this, cartMapper, isCpGrid](int elemIdx) {
300 std::array<double,dimensionworld> centroid;
301 const auto rank = this->gridView().comm().rank();
302 const auto maxLevel = this->gridView().grid().maxLevel();
303 bool useEclipse = !isCpGrid || (isCpGrid && (rank == 0) && (maxLevel == 0));
304 if (useEclipse)
305 {
306 centroid = this->eclState().getInputGrid().getCellCenter(cartMapper.cartesianIndex(elemIdx));
307 }
308 else
309 {
311 centroid = lookUpCellCentroid(elemIdx);
312 }
313 return centroid;
314 };
315 }
316
317 void callImplementationInit()
318 {
319 asImp_().createGrids_();
320 asImp_().filterConnections_();
321 std::string outputDir = Parameters::Get<Parameters::OutputDir>();
322 bool enableEclCompatFile = !Parameters::Get<Parameters::EnableOpmRstFile>();
323 asImp_().updateOutputDir_(outputDir, enableEclCompatFile);
324 const std::string& dryRunString = Parameters::Get<Parameters::EnableDryRun>();
325 asImp_().updateNOSIM_(dryRunString);
326 asImp_().finalizeInit_();
327 }
328
329 void updateCartesianToCompressedMapping_()
330 {
331 std::size_t num_cells = asImp_().grid().leafGridView().size(0);
332 is_interior_.resize(num_cells);
333
334 ElementMapper elemMapper(this->gridView(), Dune::mcmgElementLayout());
335 for (const auto& element : elements(this->gridView()))
336 {
337 const auto elemIdx = elemMapper.index(element);
338 unsigned cartesianCellIdx = cartesianIndex(elemIdx);
340 if (element.partitionType() == Dune::InteriorEntity)
341 {
342 is_interior_[elemIdx] = 1;
343 }
344 else
345 {
346 is_interior_[elemIdx] = 0;
347 }
348 }
349 }
350
351 void updateCellDepths_()
352 {
353 int numCells = this->gridView().size(/*codim=*/0);
354 cellCenterDepth_.resize(numCells);
355
356 ElementMapper elemMapper(this->gridView(), Dune::mcmgElementLayout());
357
358 const auto num_aqu_cells = this->allAquiferCells();
359
360 for(const auto& element : elements(this->gridView())) {
361 const unsigned int elemIdx = elemMapper.index(element);
362 cellCenterDepth_[elemIdx] = cellCenterDepth(element);
363
364 if (!num_aqu_cells.empty()) {
365 const unsigned int global_index = cartesianIndex(elemIdx);
366 const auto search = num_aqu_cells.find(global_index);
367 if (search != num_aqu_cells.end()) {
368 // updating the cell depth using aquifer cell depth
369 cellCenterDepth_[elemIdx] = search->second->depth;
370 }
371 }
372 }
373 }
374 void updateCellThickness_()
375 {
376 if (!this->drsdtconEnabled())
377 return;
378
379 ElementMapper elemMapper(this->gridView(), Dune::mcmgElementLayout());
380
381 int numElements = this->gridView().size(/*codim=*/0);
382 cellThickness_.resize(numElements);
383
384 for (const auto& elem : elements(this->gridView())) {
385 const unsigned int elemIdx = elemMapper.index(elem);
386 cellThickness_[elemIdx] = computeCellThickness(elem);
387 }
388 }
389
390private:
391 // computed from averaging cell corner depths
392 Scalar cellCenterDepth(const Element& element) const
393 {
394 typedef typename Element::Geometry Geometry;
395 static constexpr int zCoord = Element::dimension - 1;
396 Scalar zz = 0.0;
397
398 const Geometry& geometry = element.geometry();
399 const int corners = geometry.corners();
400 for (int i=0; i < corners; ++i)
401 zz += geometry.corner(i)[zCoord];
402
403 return zz/Scalar(corners);
404 }
405
406 Scalar computeCellThickness(const typename GridView::template Codim<0>::Entity& element) const
407 {
408 typedef typename Element::Geometry Geometry;
409 static constexpr int zCoord = Element::dimension - 1;
410 Scalar zz1 = 0.0;
411 Scalar zz2 = 0.0;
412
413 const Geometry& geometry = element.geometry();
414 // This code only works with CP-grid where the
415 // number of corners are 8 and
416 // also assumes that the first
417 // 4 corners are the top surface and
418 // the 4 next are the bottomn.
419 assert(geometry.corners() == 8);
420 for (int i=0; i < 4; ++i){
421 zz1 += geometry.corner(i)[zCoord];
422 zz2 += geometry.corner(i+4)[zCoord];
423 }
424 zz1 /=4;
425 zz2 /=4;
426 return zz2-zz1;
427 }
428
429 Implementation& asImp_()
430 { return *static_cast<Implementation*>(this); }
431
432 const Implementation& asImp_() const
433 { return *static_cast<const Implementation*>(this); }
434
435protected:
439 std::unordered_map<int,int> cartesianToCompressed_;
440
443 std::vector<Scalar> cellCenterDepth_;
444
447 std::vector<Scalar> cellThickness_;
448
451 std::vector<int> is_interior_;
452};
453
454} // namespace Opm
455
456#endif // OPM_FLOW_BASE_VANGUARD_HPP
Helper class for grid instantiation of ECL file-format using problems.
Provides the base class for most (all?) simulator vanguards.
Definition CollectDataOnIORank.hpp:49
Provides the base class for most (all?) simulator vanguards.
Definition basevanguard.hh:49
const GridView & gridView() const
Returns a reference to the grid view to be used.
Definition basevanguard.hh:69
Helper class for grid instantiation of ECL file-format using problems.
Definition FlowBaseVanguard.hpp:83
unsigned cartesianIndex(const std::array< int, dimension > &coords) const
Return the index of the cells in the logical Cartesian grid.
Definition FlowBaseVanguard.hpp:161
int cartesianSize() const
Returns the overall number of cells of the logically Cartesian grid.
Definition FlowBaseVanguard.hpp:143
int compressedIndex(int cartesianCellIdx) const
Return compressed index from cartesian index.
Definition FlowBaseVanguard.hpp:179
std::vector< Scalar > cellCenterDepth_
Cell center depths.
Definition FlowBaseVanguard.hpp:443
const std::array< int, dimension > & cartesianDimensions() const
Returns the number of logically Cartesian cells in each direction.
Definition FlowBaseVanguard.hpp:137
unsigned equilCartesianIndex(unsigned compressedEquilCellIdx) const
Returns the Cartesian cell id given an element index for the grid used for equilibration.
Definition FlowBaseVanguard.hpp:219
FlowBaseVanguard(Simulator &simulator)
Create the grid for problem data files which use the ECL file format.
Definition FlowBaseVanguard.hpp:117
int compressedIndexForInterior(int cartesianCellIdx) const
Return compressed index from cartesian index only in interior.
Definition FlowBaseVanguard.hpp:194
unsigned cartesianIndex(unsigned compressedCellIdx) const
Returns the Cartesian cell id for identifaction with ECL data.
Definition FlowBaseVanguard.hpp:155
std::vector< int > is_interior_
Whether a cells is in the interior.
Definition FlowBaseVanguard.hpp:451
void cartesianCoordinate(unsigned cellIdx, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition FlowBaseVanguard.hpp:213
std::unordered_map< int, int > cartesianToCompressed_
Mapping between cartesian and compressed cells.
Definition FlowBaseVanguard.hpp:439
std::vector< Scalar > cellThickness_
Cell thickness.
Definition FlowBaseVanguard.hpp:447
Scalar cellCenterDepth(unsigned globalSpaceIdx) const
Returns the depth of a degree of freedom [m].
Definition FlowBaseVanguard.hpp:238
std::size_t globalNumCells() const
Get the number of cells in the global leaf grid view.
Definition FlowBaseVanguard.hpp:266
void equilCartesianCoordinate(unsigned cellIdx, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell of the grid used for EQUIL.
Definition FlowBaseVanguard.hpp:228
std::function< std::array< double, dimensionworld >(int)> cellCentroids_(const CartMapper &cartMapper, const bool &isCpGrid) const
Get function to query cell centroids for a distributed grid.
Definition FlowBaseVanguard.hpp:297
Scalar cellThickness(unsigned globalSpaceIdx) const
Returns the thickness of a degree of freedom [m].
Definition FlowBaseVanguard.hpp:255
static void registerParameters()
Register the common run-time parameters for all ECL simulator vanguards.
Definition FlowBaseVanguard.hpp:106
int equilCartesianSize() const
Returns the overall number of cells of the logically EquilCartesian grid.
Definition FlowBaseVanguard.hpp:149
Definition FlowGenericVanguard.hpp:107
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:167
Declare the properties used by the infrastructure code of the finite volume discretizations.
Declare the properties used by the infrastructure code of the finite volume discretizations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Definition FlowBaseVanguard.hpp:56
Definition FlowBaseVanguard.hpp:69
Definition FlowBaseVanguard.hpp:63
a tag to mark properties as undefined
Definition propertysystem.hh:40