27#ifndef OPM_ALUGRID_VANGUARD_HPP
28#define OPM_ALUGRID_VANGUARD_HPP
30#include <dune/alugrid/common/fromtogridfactory.hh>
31#include <dune/alugrid/dgf.hh>
32#include <dune/alugrid/grid.hh>
34#include <opm/common/OpmLog/OpmLog.hpp>
36#include <opm/grid/CpGrid.hpp>
37#include <opm/grid/cpgrid/LevelCartesianIndexMapper.hpp>
41#include <opm/simulators/flow/AluGridCartesianIndexMapper.hpp>
42#include <opm/simulators/flow/AluGridLevelCartesianIndexMapper.hpp>
45#include <opm/simulators/utils/ParallelEclipseState.hpp>
54template <
class TypeTag>
59namespace Opm::Properties {
63 using InheritsFrom = std::tuple<FlowBaseVanguard>;
68template<
class TypeTag>
72template<
class TypeTag>
75 using type = Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming, Dune::ALUGridMPIComm>;
77 using type = Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming, Dune::ALUGridNoComm>;
80template<
class TypeTag>
82 using type = Dune::CpGrid;
96template <
class TypeTag>
114 using Factory = Dune::FromToGridFactory<Grid>;
116 static constexpr int dimension = Grid::dimension;
117 static constexpr int dimensionworld = Grid::dimensionworld;
123 this->callImplementationInit();
148 {
return *equilGrid_; }
159 delete equilCartesianIndexMapper_;
160 equilCartesianIndexMapper_ =
nullptr;
163 equilGrid_ =
nullptr;
174 auto dataHandle = cartesianIndexMapper_->dataHandle(
gridView);
175 grid().loadBalance(*dataHandle);
178 grid().communicate(*dataHandle,
179 Dune::InteriorBorder_All_Interface,
180 Dune::ForwardCommunication );
184 globalTrans_ = std::make_unique<TransmissibilityType>(this->
eclState(),
196 globalTrans_->update(
false, TransmissibilityType::TransUpdateQuantities::Trans,
197 [&](
unsigned int i) {
return gridEquilIdxToGridIdx(i);});
207 template<
class DataHandle>
213 template<
class DataHandle>
219 template<
class DataHandle,
class InterfaceType,
class CommunicationDirection>
233 globalTrans_.reset();
241 {
return *cartesianIndexMapper_; }
255 {
return *equilCartesianIndexMapper_; }
264 std::function<std::array<double,dimensionworld>(
int)>
270 const TransmissibilityType& globalTransmissibility()
const
272 assert( globalTrans_ !=
nullptr );
273 return *globalTrans_;
276 const std::vector<int>& globalCell()
278 return cartesianCellId_;
281 std::vector<int> cellPartition()
const
287 unsigned int gridEquilIdxToGridIdx(
unsigned int elemIndex)
const {
291 unsigned int gridIdxToEquilGridIdx(
unsigned int elemIndex)
const {
314 global_porv = this->
eclState().fieldProps().porv(
true);
315 OpmLog::info(
"\nProcessing grid");
321 this->equilGrid_ = std::make_unique<Dune::CpGrid>();
325 this->equilGrid_->processEclipseFormat(
input_grid,
331 cartesianCellId_ = this->equilGrid_->globalCell();
333 for (
unsigned i = 0; i < dimension; ++i)
334 cartesianDimension_[i] = this->equilGrid_->logicalCartesianSize()[i];
336 equilCartesianIndexMapper_ = std::make_unique<EquilCartesianIndexMapper>(*equilGrid_);
342 factory_ = std::make_unique<Factory>();
343 grid_ = factory_->convert(*equilGrid_, cartesianCellId_, ordering_);
344 OpmLog::warning(
"Space Filling Curve (SFC) ordering is enabled: see flow_blackoil_alugrid for more informations on disabling/enabling SFC reordering");
345 equilGridToGrid_.resize(ordering_.size());
346 for (std::size_t index = 0; index < ordering_.size(); ++index) {
347 equilGridToGrid_[ordering_[index]] = index;
350 cartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*grid_, cartesianDimension_, cartesianCellId_);
351 this->updateGridView_();
352 this->updateCartesianToCompressedMapping_();
353 this->updateCellDepths_();
354 this->updateCellThickness_();
357 void filterConnections_()
362 std::unique_ptr<Grid> grid_;
363 std::unique_ptr<EquilGrid> equilGrid_;
364 std::vector<int> cartesianCellId_;
365 std::vector<unsigned int> ordering_;
366 std::vector<unsigned int> equilGridToGrid_;
367 std::array<int,dimension> cartesianDimension_;
368 std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
369 std::unique_ptr<EquilCartesianIndexMapper> equilCartesianIndexMapper_;
370 std::unique_ptr<Factory> factory_;
375 std::unique_ptr<TransmissibilityType> globalTrans_;
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 AluGridVanguard.hpp:98
void loadBalance()
Distribute the simulation grid over multiple processes.
Definition AluGridVanguard.hpp:171
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition AluGridVanguard.hpp:240
const LevelCartesianIndexMapper levelCartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition AluGridVanguard.hpp:248
std::function< std::array< double, dimensionworld >(int)> cellCentroids() const
Get function to query cell centroids for a distributed grid.
Definition AluGridVanguard.hpp:265
void releaseGlobalTransmissibilities()
Free the memory occupied by the global transmissibility object.
Definition AluGridVanguard.hpp:231
const Grid & grid() const
Return a reference to the simulation grid.
Definition AluGridVanguard.hpp:135
const EquilCartesianIndexMapper & equilCartesianIndexMapper() const
Returns mapper from compressed to cartesian indices for the EQUIL grid.
Definition AluGridVanguard.hpp:254
Grid & grid()
Return a reference to the simulation grid.
Definition AluGridVanguard.hpp:129
const EquilGrid & equilGrid() const
Returns a refefence to the grid which should be used by the EQUIL initialization code.
Definition AluGridVanguard.hpp:147
void releaseEquilGrid()
Indicates that the initial condition has been computed and the memory used by the EQUIL grid can be r...
Definition AluGridVanguard.hpp:157
const GridView & gridView() const
Returns a reference to the grid view to be used.
Definition basevanguard.hh:69
Definition AluGridVanguard.hpp:55
Helper class for grid instantiation of ECL file-format using problems.
Definition FlowBaseVanguard.hpp:83
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
static Parallel::Communication & comm()
Obtain global communicator.
Definition FlowGenericVanguard.hpp:332
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:167
Definition RelpermDiagnostics.hpp:31
Definition Transmissibility.hpp:54
Defines the common properties required by the porous medium multi-phase models.
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
Enable diffusive fluxes?
Definition multiphasebaseproperties.hh:79
Enable dispersive fluxes?
Definition multiphasebaseproperties.hh:82
Specify whether energy should be considered as a conservation quantity or not.
Definition multiphasebaseproperties.hh:76
Definition FlowBaseVanguard.hpp:69
The type of the DUNE grid.
Definition basicproperties.hh:100
Definition AluGridVanguard.hpp:62
Property which provides a Vanguard (manages grids)
Definition basicproperties.hh:96