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>
40#include <opm/simulators/flow/AluGridCartesianIndexMapper.hpp>
43#include <opm/simulators/utils/ParallelEclipseState.hpp>
52template <
class TypeTag>
57namespace Opm::Properties {
61 using InheritsFrom = std::tuple<FlowBaseVanguard>;
66template<
class TypeTag>
70template<
class TypeTag>
73 using type = Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming, Dune::ALUGridMPIComm>;
75 using type = Dune::ALUGrid<3, 3, Dune::cube, Dune::nonconforming, Dune::ALUGridNoComm>;
78template<
class TypeTag>
80 using type = Dune::CpGrid;
94template <
class TypeTag>
111 using Factory = Dune::FromToGridFactory<Grid>;
113 static constexpr int dimension = Grid::dimension;
114 static constexpr int dimensionworld = Grid::dimensionworld;
120 this->callImplementationInit();
145 {
return *equilGrid_; }
156 delete equilCartesianIndexMapper_;
157 equilCartesianIndexMapper_ =
nullptr;
160 equilGrid_ =
nullptr;
171 auto dataHandle = cartesianIndexMapper_->dataHandle(
gridView);
172 grid().loadBalance(*dataHandle);
175 grid().communicate(*dataHandle,
176 Dune::InteriorBorder_All_Interface,
177 Dune::ForwardCommunication );
181 globalTrans_ = std::make_unique<TransmissibilityType>(this->
eclState(),
193 globalTrans_->update(
false, TransmissibilityType::TransUpdateQuantities::Trans,
194 [&](
unsigned int i) {
return gridEquilIdxToGridIdx(i);});
199 template<
class DataHandle>
200 void scatterData(DataHandle& )
const
205 template<
class DataHandle>
206 void gatherData(DataHandle& )
const
211 template<
class DataHandle,
class InterfaceType,
class CommunicationDirection>
212 void communicate (DataHandle& , InterfaceType ,
213 CommunicationDirection )
const
225 globalTrans_.reset();
233 {
return *cartesianIndexMapper_; }
239 {
return *equilCartesianIndexMapper_; }
248 std::function<std::array<double,dimensionworld>(
int)>
254 const TransmissibilityType& globalTransmissibility()
const
256 assert( globalTrans_ !=
nullptr );
257 return *globalTrans_;
260 const std::vector<int>& globalCell()
262 return cartesianCellId_;
265 std::vector<int> cellPartition()
const
271 unsigned int gridEquilIdxToGridIdx(
unsigned int elemIndex)
const {
272 return equilGridToGrid_[elemIndex];
275 unsigned int gridIdxToEquilGridIdx(
unsigned int elemIndex)
const {
276 return ordering_[elemIndex];
290 const EclipseGrid* input_grid =
nullptr;
291 std::vector<double> global_porv;
297 input_grid = &this->
eclState().getInputGrid();
298 global_porv = this->
eclState().fieldProps().porv(
true);
299 OpmLog::info(
"\nProcessing grid");
305 this->equilGrid_ = std::make_unique<Dune::CpGrid>();
309 this->equilGrid_->processEclipseFormat(input_grid,
315 cartesianCellId_ = this->equilGrid_->globalCell();
317 for (
unsigned i = 0; i < dimension; ++i)
318 cartesianDimension_[i] = this->equilGrid_->logicalCartesianSize()[i];
320 equilCartesianIndexMapper_ = std::make_unique<EquilCartesianIndexMapper>(*equilGrid_);
326 factory_ = std::make_unique<Factory>();
327 grid_ = factory_->convert(*equilGrid_, cartesianCellId_, ordering_);
328 OpmLog::warning(
"Space Filling Curve (SFC) ordering is enabled: see flow_blackoil_alugrid for more informations on disabling/enabling SFC reordering");
329 equilGridToGrid_.resize(ordering_.size());
330 for (std::size_t index = 0; index < ordering_.size(); ++index) {
331 equilGridToGrid_[ordering_[index]] = index;
334 cartesianIndexMapper_ = std::make_unique<CartesianIndexMapper>(*grid_, cartesianDimension_, cartesianCellId_);
335 this->updateGridView_();
336 this->updateCartesianToCompressedMapping_();
337 this->updateCellDepths_();
338 this->updateCellThickness_();
341 void filterConnections_()
346 std::unique_ptr<Grid> grid_;
347 std::unique_ptr<EquilGrid> equilGrid_;
348 std::vector<int> cartesianCellId_;
349 std::vector<unsigned int> ordering_;
350 std::vector<unsigned int> equilGridToGrid_;
351 std::array<int,dimension> cartesianDimension_;
352 std::unique_ptr<CartesianIndexMapper> cartesianIndexMapper_;
353 std::unique_ptr<EquilCartesianIndexMapper> equilCartesianIndexMapper_;
354 std::unique_ptr<Factory> factory_;
359 std::unique_ptr<TransmissibilityType> globalTrans_;
Helper class for grid instantiation of ECL file-format using problems.
Definition findOverlapRowsAndColumns.hpp:31
Helper class for grid instantiation of ECL file-format using problems.
Definition AluGridVanguard.hpp:96
void loadBalance()
Distribute the simulation grid over multiple processes.
Definition AluGridVanguard.hpp:168
const CartesianIndexMapper & cartesianIndexMapper() const
Returns the object which maps a global element index of the simulation grid to the corresponding elem...
Definition AluGridVanguard.hpp:232
std::function< std::array< double, dimensionworld >(int)> cellCentroids() const
Get function to query cell centroids for a distributed grid.
Definition AluGridVanguard.hpp:249
void releaseGlobalTransmissibilities()
Free the memory occupied by the global transmissibility object.
Definition AluGridVanguard.hpp:223
const Grid & grid() const
Return a reference to the simulation grid.
Definition AluGridVanguard.hpp:132
const EquilCartesianIndexMapper & equilCartesianIndexMapper() const
Returns mapper from compressed to cartesian indices for the EQUIL grid.
Definition AluGridVanguard.hpp:238
Grid & grid()
Return a reference to the simulation grid.
Definition AluGridVanguard.hpp:126
const EquilGrid & equilGrid() const
Returns a refefence to the grid which should be used by the EQUIL initialization code.
Definition AluGridVanguard.hpp:144
void releaseEquilGrid()
Indicates that the initial condition has been computed and the memory used by the EQUIL grid can be r...
Definition AluGridVanguard.hpp:154
const GridView & gridView() const
Returns a reference to the grid view to be used.
Definition basevanguard.hh:69
Definition AluGridVanguard.hpp:53
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:318
static Parallel::Communication & comm()
Obtain global communicator.
Definition FlowGenericVanguard.hpp:306
const EclipseState & eclState() const
Return a reference to the internalized ECL deck.
Definition FlowGenericVanguard.hpp:159
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
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
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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:60
Property which provides a Vanguard (manages grids)
Definition basicproperties.hh:96