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>();
120 fileName_ = Parameters::Get<Parameters::EclDeckFileName>();
121 edgeWeightsMethod_ = Dune::EdgeWeightMethod(Parameters::Get<Parameters::EdgeWeightsMethod>());
123#if HAVE_OPENCL || HAVE_ROCSPARSE || HAVE_CUDA
124 numJacobiBlocks_ = Parameters::Get<Parameters::NumJacobiBlocks>();
127 ownersFirst_ = Parameters::Get<Parameters::OwnerCellsFirst>();
129 partitionMethod_ = Dune::PartitionMethod(Parameters::Get<Parameters::PartitionMethod>());
131 imbalanceTol_ = Parameters::Get<Parameters::ImbalanceTol<Scalar>>();
137 metisParams_ = Parameters::Get<Parameters::MetisParams>();
141 enableDistributedWells_ = Parameters::Get<Parameters::AllowDistributedWells>();
142 ignoredKeywords_ = Parameters::Get<Parameters::IgnoreKeywords>();
143 int output_param = Parameters::Get<Parameters::EclOutputInterval>();
146 useMultisegmentWell_ = Parameters::Get<Parameters::UseMultisegmentWell>();
147 enableExperiments_ = enableExperiments;
152 const CartesianIndexMapper& cartesianMapper()
const
153 {
return asImp_().cartesianIndexMapper(); }
159 {
return asImp_().cartesianIndexMapper().cartesianDimensions(); }
165 {
return asImp_().cartesianIndexMapper().cartesianSize(); }
171 {
return asImp_().equilCartesianIndexMapper().cartesianSize(); }
186 for (
unsigned i = 1; i < dimension; ++i) {
235 {
return asImp_().cartesianIndexMapper().cartesianCoordinate(
cellIdx,
ijk); }
250 {
return asImp_().equilCartesianIndexMapper().cartesianCoordinate(
cellIdx,
ijk); }
264 const std::vector<Scalar>& cellCenterDepths()
const
289 const auto& grid = asImp_().grid();
290 if (grid.comm().size() == 1)
292 return grid.leafGridView().size(0);
294 const auto&
gridView = grid.leafGridView();
295 constexpr int codim = 0;
296 constexpr auto Part = Dune::Interior_Partition;
302 void setupCartesianToCompressed_() {
303 this->updateCartesianToCompressedMapping_();
316 template<
class CartMapper>
317 std::function<std::array<double,dimensionworld>(
int)>
321 std::array<double,dimensionworld>
centroid;
322 const auto rank = this->
gridView().comm().rank();
338 void callImplementationInit()
340 asImp_().createGrids_();
341 asImp_().filterConnections_();
342 std::string outputDir = Parameters::Get<Parameters::OutputDir>();
345 const std::string&
dryRunString = Parameters::Get<Parameters::EnableDryRun>();
347 asImp_().finalizeInit_();
350 void updateCartesianToCompressedMapping_()
352 std::size_t
num_cells = asImp_().grid().leafGridView().size(0);
361 if (element.partitionType() == Dune::InteriorEntity)
372 void updateCellDepths_()
374 int numCells = this->
gridView().size(0);
395 void updateCellThickness_()
397 if (!this->drsdtconEnabled())
415 typedef typename Element::Geometry Geometry;
416 static constexpr int zCoord = Element::dimension - 1;
419 const Geometry& geometry = element.geometry();
420 const int corners = geometry.corners();
421 for (
int i=0; i <
corners; ++i)
427 Scalar computeCellThickness(
const typename GridView::template Codim<0>::Entity& element)
const
429 typedef typename Element::Geometry Geometry;
430 static constexpr int zCoord = Element::dimension - 1;
434 const Geometry& geometry = element.geometry();
440 assert(geometry.corners() == 8);
441 for (
int i=0; i < 4; ++i){
450 Implementation& asImp_()
451 {
return *
static_cast<Implementation*
>(
this); }
453 const Implementation& asImp_()
const
454 {
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:51
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:182
int cartesianSize() const
Returns the overall number of cells of the logically Cartesian grid.
Definition FlowBaseVanguard.hpp:164
int compressedIndex(int cartesianCellIdx) const
Return compressed index from cartesian index.
Definition FlowBaseVanguard.hpp:200
std::vector< Scalar > cellCenterDepth_
Cell center depths.
Definition FlowBaseVanguard.hpp:464
const std::array< int, dimension > & cartesianDimensions() const
Returns the number of logically Cartesian cells in each direction.
Definition FlowBaseVanguard.hpp:158
unsigned equilCartesianIndex(unsigned compressedEquilCellIdx) const
Returns the Cartesian cell id given an element index for the grid used for equilibration.
Definition FlowBaseVanguard.hpp:240
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:215
unsigned cartesianIndex(unsigned compressedCellIdx) const
Returns the Cartesian cell id for identifaction with ECL data.
Definition FlowBaseVanguard.hpp:176
std::vector< int > is_interior_
Whether a cells is in the interior.
Definition FlowBaseVanguard.hpp:472
void cartesianCoordinate(unsigned cellIdx, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition FlowBaseVanguard.hpp:234
std::unordered_map< int, int > cartesianToCompressed_
Mapping between cartesian and compressed cells.
Definition FlowBaseVanguard.hpp:460
std::vector< Scalar > cellThickness_
Cell thickness.
Definition FlowBaseVanguard.hpp:468
Scalar cellCenterDepth(unsigned globalSpaceIdx) const
Returns the depth of a degree of freedom [m].
Definition FlowBaseVanguard.hpp:259
std::size_t globalNumCells() const
Get the number of cells in the global leaf grid view.
Definition FlowBaseVanguard.hpp:287
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:249
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:318
Scalar cellThickness(unsigned globalSpaceIdx) const
Returns the thickness of a degree of freedom [m].
Definition FlowBaseVanguard.hpp:276
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:170
Definition FlowGenericVanguard.hpp:99
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:159
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