27#ifndef OPM_POLYHEDRAL_GRID_VANGUARD_HPP
28#define OPM_POLYHEDRAL_GRID_VANGUARD_HPP
30#include <opm/grid/polyhedralgrid.hh>
32#include <opm/models/common/multiphasebaseproperties.hh>
41#include <unordered_set>
44template <
class TypeTag>
45class PolyhedralGridVanguard;
48namespace Opm::Properties {
52 using InheritsFrom = std::tuple<FlowBaseVanguard>;
57template<
class TypeTag>
61template<
class TypeTag>
63 using type = Dune::PolyhedralGrid<3, 3>;
65template<
class TypeTag>
67 using type = GetPropType<TypeTag, Properties::Grid>;
81template <
class TypeTag>
87 using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
88 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
89 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
92 using Grid = GetPropType<TypeTag, Properties::Grid>;
93 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
94 using GridView = GetPropType<TypeTag, Properties::GridView>;
97 static constexpr int dimension = Grid::dimension;
98 static constexpr int dimensionworld = Grid::dimensionworld;
101 using GridPointer = Grid*;
102 using EquilGridPointer = EquilGrid*;
105 using TransmissibilityType =
Transmissibility<Grid, GridView, ElementMapper,
110 , simulator_(simulator)
112 this->callImplementationInit();
114 const int* globalcellorg = this->
grid().globalCell();
115 int num_cells = this->gridView().size(0);
116 globalcell_.resize(num_cells);
117 for(
int i=0; i < num_cells; ++i){
118 globalcell_[i] = globalcellorg[i];
170 {
return *cartesianIndexMapper_; }
179 {
return *cartesianIndexMapper_; }
181 const std::vector<int>& globalCell()
186 unsigned int gridEquilIdxToGridIdx(
unsigned int elemIndex)
const {
190 unsigned int gridIdxToEquilGridIdx(
unsigned int elemIndex)
const {
203 std::unordered_set<std::string> defunctWellNames()
const
204 {
return defunctWellNames_; }
206 const TransmissibilityType& globalTransmissibility()
const
208 return simulator_.problem().eclTransmissibilities();
218 std::function<std::array<double,FlowBaseVanguard<TypeTag>::dimensionworld>(int)>
224 std::vector<int> cellPartition()
const
232 grid_ = std::make_unique<Grid>(this->
eclState().getInputGrid(), this->
eclState().fieldProps().porv(
true));
233 cartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*grid_);
234 this->updateGridView_();
235 this->updateCartesianToCompressedMapping_();
236 this->updateCellDepths_();
239 void filterConnections_()
244 Simulator& simulator_;
246 std::unique_ptr<Grid> grid_;
247 std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
250 std::unordered_set<std::string> defunctWellNames_;
251 std::vector<int> globalcell_;
Helper class for grid instantiation of ECL file-format using problems.
Definition CollectDataOnIORank.hpp:49
Helper class for grid instantiation of ECL file-format using problems.
Definition FlowBaseVanguard.hpp:216
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:496
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:120
Helper class for grid instantiation of ECL file-format using problems.
Definition PolyhedralGridVanguard.hpp:83
void releaseEquilGrid()
Indicates that the initial condition has been computed and the memory used by the EQUIL grid can be r...
Definition PolyhedralGridVanguard.hpp:153
const CartesianIndexMapper & equilCartesianIndexMapper() const
Returns mapper from compressed to cartesian indices for the EQUIL grid.
Definition PolyhedralGridVanguard.hpp:178
void releaseGlobalTransmissibilities()
Free the memory occupied by the global transmissibility object.
Definition PolyhedralGridVanguard.hpp:199
void loadBalance()
Distribute the simulation grid over multiple processes.
Definition PolyhedralGridVanguard.hpp:161
const EquilGrid & equilGrid() const
Returns a refefence to the grid which should be used by the EQUIL initialization code.
Definition PolyhedralGridVanguard.hpp:143
Grid & grid()
Return a reference to the simulation grid.
Definition PolyhedralGridVanguard.hpp:125
const Grid & grid() const
Return a reference to the simulation grid.
Definition PolyhedralGridVanguard.hpp:131
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition PolyhedralGridVanguard.hpp:169
std::function< std::array< double, FlowBaseVanguard< TypeTag >::dimensionworld >(int)> cellCentroids() const
Get function to query cell centroids for a distributed grid.
Definition PolyhedralGridVanguard.hpp:219
Definition Transmissibility.hpp:54
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition FlowBaseVanguard.hpp:64
Definition PolyhedralGridVanguard.hpp:51