27#ifndef OPM_FLOW_BASE_VANGUARD_HPP
28#define OPM_FLOW_BASE_VANGUARD_HPP
30#include <dune/grid/common/partitionset.hh>
32#include <opm/grid/common/GridEnums.hpp>
33#include <opm/grid/common/CartesianIndexMapper.hpp>
34#include <opm/grid/LookUpCellCentroid.hh>
36#include <opm/input/eclipse/EclipseState/Aquifer/NumericalAquifer/NumericalAquiferCell.hpp>
37#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
45#include <opm/simulators/flow/BlackoilModelParameters.hpp>
50#include <unordered_map>
54template <
class TypeTag>
55class FlowBaseVanguard;
59namespace Opm::Properties {
68template<
class TypeTag,
class MyTypeTag>
80template <
class TypeTag>
97 static const int dimension = Grid::dimension;
98 static const int dimensionworld = Grid::dimensionworld;
99 using Element =
typename GridView::template Codim<0>::Entity;
108 FlowGenericVanguard::registerParameters_<Scalar>();
121 imbalanceTol_ = Parameters::Get<Parameters::ImbalanceTol<Scalar>>();
126 enableExperiments_ = enableExperiments;
131 const CartesianIndexMapper& cartesianMapper()
const
132 {
return asImp_().cartesianIndexMapper(); }
138 {
return asImp_().cartesianIndexMapper().cartesianDimensions(); }
144 {
return asImp_().cartesianIndexMapper().cartesianSize(); }
150 {
return asImp_().equilCartesianIndexMapper().cartesianSize(); }
165 for (
unsigned i = 1; i < dimension; ++i) {
214 {
return asImp_().cartesianIndexMapper().cartesianCoordinate(
cellIdx,
ijk); }
229 {
return asImp_().equilCartesianIndexMapper().cartesianCoordinate(
cellIdx,
ijk); }
243 const std::vector<Scalar>& cellCenterDepths()
const
268 const auto& grid = asImp_().grid();
269 if (grid.comm().size() == 1)
271 return grid.leafGridView().size(0);
273 const auto&
gridView = grid.leafGridView();
274 constexpr int codim = 0;
275 constexpr auto Part = Dune::Interior_Partition;
281 void setupCartesianToCompressed_() {
282 this->updateCartesianToCompressedMapping_();
295 template<
class CartMapper>
296 std::function<std::array<double,dimensionworld>(
int)>
300 std::array<double,dimensionworld>
centroid;
301 const auto rank = this->
gridView().comm().rank();
306 centroid = this->
eclState().getInputGrid().getCellCenter(cartMapper.cartesianIndex(elemIdx));
317 void callImplementationInit()
319 asImp_().createGrids_();
320 asImp_().filterConnections_();
321 std::string outputDir = Parameters::Get<Parameters::OutputDir>();
324 const std::string&
dryRunString = Parameters::Get<Parameters::EnableDryRun>();
326 asImp_().finalizeInit_();
329 void updateCartesianToCompressedMapping_()
331 std::size_t
num_cells = asImp_().grid().leafGridView().size(0);
337 const auto elemIdx =
elemMapper.index(element);
340 if (element.partitionType() == Dune::InteriorEntity)
351 void updateCellDepths_()
353 int numCells = this->
gridView().size(0);
361 const unsigned int elemIdx =
elemMapper.index(element);
374 void updateCellThickness_()
376 if (!this->drsdtconEnabled())
394 typedef typename Element::Geometry Geometry;
395 static constexpr int zCoord = Element::dimension - 1;
398 const Geometry& geometry = element.geometry();
399 const int corners = geometry.corners();
400 for (
int i=0; i <
corners; ++i)
406 Scalar computeCellThickness(
const typename GridView::template Codim<0>::Entity& element)
const
408 typedef typename Element::Geometry Geometry;
409 static constexpr int zCoord = Element::dimension - 1;
413 const Geometry& geometry = element.geometry();
419 assert(geometry.corners() == 8);
420 for (
int i=0; i < 4; ++i){
429 Implementation& asImp_()
430 {
return *
static_cast<Implementation*
>(
this); }
432 const Implementation& asImp_()
const
433 {
return *
static_cast<const Implementation*
>(
this); }
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