[ provides Dune::Grid ] More...

#include <CpGrid.hpp>

Inheritance diagram for Dune::CpGrid:
Inheritance graph

Classes

class  CentroidIterator
 An iterator over the centroids of the geometry of the entities. More...
 

Public Types

typedef CpGridFamily GridFamily
 Family typedef, why is this not defined by Grid<>? More...
 

Public Member Functions

 CpGrid ()
 Default constructor. More...
 
 CpGrid (MPIHelper::MPICommunicator comm)
 
IO routines
void readSintefLegacyFormat (const std::string &grid_prefix)
 
void writeSintefLegacyFormat (const std::string &grid_prefix) const
 
void processEclipseFormat (const grdecl &input_data, bool remove_ij_boundary, bool turn_normals=false)
 
Cartesian grid extensions.

A cornerpoint grid can be seen as a degenerated and distorted cartesian grid. Therefore it provides mappings from cells to the underlying cartesian index.

void createCartesian (const std::array< int, 3 > &dims, const std::array< double, 3 > &cellsize, const std::array< int, 3 > &shift={0, 0, 0})
 
const std::array< int, 3 > & logicalCartesianSize () const
 
const std::vector< int > & globalCell () const
 
const std::vector< std::shared_ptr< Dune::cpgrid::CpGridData > > & chooseData () const
 Returns either data_ or distributed_data_(if non empty). More...
 
void getIJK (const int c, std::array< int, 3 > &ijk) const
 Extract Cartesian index triplet (i,j,k) of an active cell. More...
 
bool uniqueBoundaryIds () const
 
void setUniqueBoundaryIds (bool uids)
 
The DUNE grid interface implementation
std::string name () const
 Get the grid name. More...
 
int maxLevel () const
 Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview), 0 = the coarsest level. More...
 
template<int codim>
Traits::template Codim< codim >::LevelIterator lbegin (int level) const
 Iterator to first entity of given codim on level. More...
 
template<int codim>
Traits::template Codim< codim >::LevelIterator lend (int level) const
 one past the end on this level More...
 
template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin (int level) const
 Iterator to first entity of given codim on level and PartitionIteratorType. More...
 
template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lend (int level) const
 one past the end on this level More...
 
template<int codim>
Traits::template Codim< codim >::LeafIterator leafbegin () const
 Iterator to first leaf entity of given codim. More...
 
template<int codim>
Traits::template Codim< codim >::LeafIterator leafend () const
 one past the end of the sequence of leaf entities More...
 
template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin () const
 Iterator to first leaf entity of given codim and PartitionIteratorType. More...
 
template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend () const
 one past the end of the sequence of leaf entities More...
 
int size (int level, int codim) const
 Number of grid entities per level and codim. More...
 
int size (int codim) const
 number of leaf entities per codim in this process More...
 
int size (int level, GeometryType type) const
 number of entities per level and geometry type in this process More...
 
int size (GeometryType type) const
 number of leaf entities per geometry type in this process More...
 
const Traits::GlobalIdSet & globalIdSet () const
 Access to the GlobalIdSet. More...
 
const Traits::LocalIdSet & localIdSet () const
 Access to the LocalIdSet. More...
 
const Traits::LevelIndexSet & levelIndexSet (int level) const
 Access to the LevelIndexSets. More...
 
const Traits::LeafIndexSet & leafIndexSet () const
 Access to the LeafIndexSet. More...
 
void globalRefine (int)
 global refinement More...
 
const std::vector< Dune::GeometryType > & geomTypes (const int) const
 
template<int codim>
cpgrid::Entity< codim > entity (const cpgrid::Entity< codim > &seed) const
 given an EntitySeed (or EntityPointer) return an entity object More...
 
void addLgrUpdateLeafView (const std::array< int, 3 > &cells_per_dim, const std::array< int, 3 > &startIJK, const std::array< int, 3 > &endIJK, const std::string &lgr_name)
 Create a grid out of a coarse one and a refinement(LGR) of a selected block-shaped patch of cells from that coarse grid. More...
 
void addLgrsUpdateLeafView (const std::vector< std::array< int, 3 > > &cells_per_dim_vec, const std::vector< std::array< int, 3 > > &startIJK_vec, const std::vector< std::array< int, 3 > > &endIJK_vec, const std::vector< std::string > &lgr_name_vec)
 Create a grid out of a coarse one and (at most) 2 refinements(LGRs) of selected block-shaped disjoint patches of cells from that coarse grid. More...
 
const std::map< std::string, int > & getLgrNameToLevel () const
 
std::array< double, 3 > getEclCentroid (const int &idx) const
 
std::array< double, 3 > getEclCentroid (const cpgrid::Entity< 0 > &elem) const
 
Dune::cpgrid::Intersection getParentIntersectionFromLgrBoundaryFace (const Dune::cpgrid::Intersection &intersection) const
 
unsigned int overlapSize (int) const
 Size of the overlap on the leaf level. More...
 
unsigned int ghostSize (int) const
 Size of the ghost cell layer on the leaf level. More...
 
unsigned int overlapSize (int, int) const
 Size of the overlap on a given level. More...
 
unsigned int ghostSize (int, int) const
 Size of the ghost cell layer on a given level. More...
 
unsigned int numBoundarySegments () const
 returns the number of boundary segments within the macro grid More...
 
void setZoltanParams (const std::map< std::string, std::string > &params)
 
bool loadBalance (int overlapLayers=1, bool useZoltan=true)
 Distributes this grid over the available nodes in a distributed machine. More...
 
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance (const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, int overlapLayers=1, bool useZoltan=true)
 Distributes this grid over the available nodes in a distributed machine. More...
 
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance (EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, bool useZoltan=true)
 Distributes this grid over the available nodes in a distributed machine. More...
 
template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance (DataHandle &data, const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities=nullptr, int overlapLayers=1, bool useZoltan=true)
 Distributes this grid and data over the available nodes in a distributed machine. More...
 
template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance (DataHandle &data, EdgeWeightMethod method, const std::vector< cpgrid::OpmWellType > *wells, bool serialPartitioning, const double *transmissibilities=nullptr, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1, bool useZoltan=true, double zoltanImbalanceTol=1.1, bool allowDistributedWells=false)
 Distributes this grid over the available nodes in a distributed machine. More...
 
template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > loadBalance (DataHandle &data, const std::vector< int > &parts, const std::vector< cpgrid::OpmWellType > *wells, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
 Distributes this grid over the available nodes in a distributed machine. More...
 
template<class DataHandle >
bool loadBalance (DataHandle &data, decltype(data.fixedSize(0, 0)) overlapLayers=1, bool useZoltan=true)
 Distributes this grid and data over the available nodes in a distributed machine. More...
 
bool loadBalance (const std::vector< int > &parts, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
 Distributes this grid over the available nodes in a distributed machine. More...
 
template<class DataHandle >
bool loadBalance (DataHandle &data, const std::vector< int > &parts, bool ownersFirst=false, bool addCornerCells=false, int overlapLayers=1)
 Distributes this grid and data over the available nodes in a distributed machine. More...
 
std::vector< int > zoltanPartitionWithoutScatter (const std::vector< cpgrid::OpmWellType > *wells, const double *transmissibilities, const int numParts, const double zoltanImbalanceTol) const
 Partitions the grid using Zoltan without decomposing and distributing it among processes. More...
 
template<class DataHandle >
void communicate (DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
 communicate objects for all codims on a given level More...
 
template<class DataHandle >
void communicate (DataHandle &data, InterfaceType iftype, CommunicationDirection dir) const
 communicate objects for all codims on a given level. More...
 
const CpGridTraits::Communicationcomm () const
 Get the collective communication object. More...
 

Friends

class cpgrid::CpGridData
 
class cpgrid::Entity< 0 >
 
class cpgrid::Entity< 1 >
 
class cpgrid::Entity< 2 >
 
class cpgrid::Entity< 3 >
 
template<typename Grid , typename GridView >
class Opm::LookUpData
 
template<typename Grid , typename GridView >
class Opm::LookUpCartesianData
 
template<int dim>
cpgrid::Entity< dim > createEntity (const CpGrid &, int, bool)
 

The simplified grid interface.

It provides additional methods not in the DUNE interface but needed by OPM. Vertices, faces, and cells are not represented as entities but identified by indices.

typedef Dune::FieldVector< double, 3 > Vector
 
const std::vector< double > & zcornData () const
 
int numCells () const
 Get the number of cells. More...
 
int numFaces () const
 Get the number of faces. More...
 
int numVertices () const
 Get The number of vertices. More...
 
int numCellFaces (int cell) const
 Get the number of faces of a cell. More...
 
int cellFace (int cell, int local_index) const
 Get a specific face of a cell. More...
 
const cpgrid::OrientedEntityTable< 0, 1 >::row_type cellFaceRow (int cell) const
 Get a list of indices identifying all faces of a cell. More...
 
int faceCell (int face, int local_index) const
 Get the index identifying a cell attached to a face. More...
 
int numCellFaces () const
 Get the sum of all faces attached to all cells. More...
 
int numFaceVertices (int face) const
 
int faceVertex (int face, int local_index) const
 Get the index identifying a vertex of a face. More...
 
double cellCenterDepth (int cell_index) const
 Get vertical position of cell center ("zcorn" average). More...
 
const Vector faceCenterEcl (int cell_index, int face, const Dune::cpgrid::Intersection &intersection) const
 
const Vector faceAreaNormalEcl (int face) const
 
const VectorvertexPosition (int vertex) const
 Get the Position of a vertex. More...
 
double faceArea (int face) const
 Get the area of a face. More...
 
const VectorfaceCentroid (int face) const
 Get the coordinates of the center of a face. More...
 
const VectorfaceNormal (int face) const
 Get the unit normal of a face. More...
 
double cellVolume (int cell) const
 Get the volume of the cell. More...
 
const VectorcellCentroid (int cell) const
 Get the coordinates of the center of a cell. More...
 
CentroidIterator< 0 > beginCellCentroids () const
 Get an iterator over the cell centroids positioned at the first one. More...
 
CentroidIterator< 1 > beginFaceCentroids () const
 Get an iterator over the face centroids positioned at the first one. More...
 
int boundaryId (int face) const
 
template<class Cell2FacesRowIterator >
int faceTag (const Cell2FacesRowIterator &cell_face) const
 Get the cartesian tag associated with a face tag. More...
 

Parallel grid extensions.

Methods extending the DUNE's parallel grid interface. These are basically for scattering/gathering data to/from distributed views.

using InterfaceMap = cpgrid::CpGridDataTraits::InterfaceMap
 The type of the map describing communication interfaces. More...
 
using ParallelIndexSet = cpgrid::CpGridDataTraits::ParallelIndexSet
 The type of the parallel index set. More...
 
using RemoteIndices = cpgrid::CpGridDataTraits::RemoteIndices
 The type of the remote indices information. More...
 
using CommunicationType = cpgrid::CpGridDataTraits::CommunicationType
 The type of the owner-overlap-copy communication. More...
 
template<class DataHandle >
void scatterData (DataHandle &handle) const
 Moves data from the global (all data on process) view to the distributed view. More...
 
template<class DataHandle >
void gatherData (DataHandle &handle) const
 Moves data from the distributed view to the global (all data on process) view. More...
 
const InterfaceMapcellScatterGatherInterface () const
 Get an interface for gathering/scattering data attached to cells with communication. More...
 
const InterfaceMappointScatterGatherInterface () const
 Get an interface for gathering/scattering data attached to points with communication. More...
 
void switchToGlobalView ()
 Switch to the global view. More...
 
void switchToDistributedView ()
 Switch to the distributed view. More...
 
const CommunicationTypecellCommunication () const
 Get the owner-overlap-copy communication for cells. More...
 
ParallelIndexSetgetCellIndexSet ()
 
RemoteIndicesgetCellRemoteIndices ()
 
const ParallelIndexSetgetCellIndexSet () const
 
const RemoteIndicesgetCellRemoteIndices () const
 
const std::vector< int > & sortedNumAquiferCells () const
 Get sorted active cell indices of numerical aquifer. More...
 

Detailed Description

[ provides Dune::Grid ]

Member Typedef Documentation

◆ CommunicationType

The type of the owner-overlap-copy communication.

◆ GridFamily

Family typedef, why is this not defined by Grid<>?

◆ InterfaceMap

The type of the map describing communication interfaces.

◆ ParallelIndexSet

◆ RemoteIndices

The type of the remote indices information.

◆ Vector

typedef Dune::FieldVector<double, 3> Dune::CpGrid::Vector

Constructor & Destructor Documentation

◆ CpGrid() [1/2]

Dune::CpGrid::CpGrid ( )

Default constructor.

◆ CpGrid() [2/2]

Dune::CpGrid::CpGrid ( MPIHelper::MPICommunicator  comm)

Member Function Documentation

◆ addLgrsUpdateLeafView()

void Dune::CpGrid::addLgrsUpdateLeafView ( const std::vector< std::array< int, 3 > > &  cells_per_dim_vec,
const std::vector< std::array< int, 3 > > &  startIJK_vec,
const std::vector< std::array< int, 3 > > &  endIJK_vec,
const std::vector< std::string > &  lgr_name_vec 
)

Create a grid out of a coarse one and (at most) 2 refinements(LGRs) of selected block-shaped disjoint patches of cells from that coarse grid.

Level0 refers to the coarse grid, assumed to be this-> data_[0]. Level1 and level2 refer to the LGRs (stored in this->data_[1] data_[2]). LeafView (stored in this-> data_[3]) is built with the level0-entities which weren't involded in the refinenment, together with the new born entities created in level1 and level2. Old-corners and old-faces (from coarse grid) lying on the boundary of the patches, get replaced by new-born-equivalent corners and new-born-faces.

Parameters
[in]cells_per_dim_vecVector of Number of (refined) cells in each direction that each parent cell should be refined to.
[in]startIJK_vecVector of Cartesian triplet indices where each patch starts.
[in]endIJK_vecVector of Cartesian triplet indices where each patch ends. Last cell part of each patch(lgr) will be {endIJK_vec[<patch-number>][0]-1, ..., endIJK_vec[<patch-number>][2]-1}.
[in]lgr_name_vecNames (std::string) for the LGRs/levels

◆ addLgrUpdateLeafView()

void Dune::CpGrid::addLgrUpdateLeafView ( const std::array< int, 3 > &  cells_per_dim,
const std::array< int, 3 > &  startIJK,
const std::array< int, 3 > &  endIJK,
const std::string &  lgr_name 
)

Create a grid out of a coarse one and a refinement(LGR) of a selected block-shaped patch of cells from that coarse grid.

Level0 refers to the coarse grid, assumed to be this-> data_[0]. Level1 refers to the LGR (stored in this->data_[1]). LeafView (stored in this-> data_[2]) is built with the level0-entities which weren't involded in the refinenment, together with the new born entities created in level1. Old-corners and old-faces (from coarse grid) lying on the boundary of the patch, get replaced by new-born-equivalent corners and new-born-faces.

Parameters
[in]cells_per_dimNumber of (refined) cells in each direction that each parent cell should be refined to.
[in]startIJKCartesian triplet index where the patch starts.
[in]endIJKCartesian triplet index where the patch ends. Last cell part of the lgr will be {endijk[0]-1, ... endIJK[2]-1}.
[in]lgr_nameName (std::string) for the lgr/level1

◆ beginCellCentroids()

CentroidIterator< 0 > Dune::CpGrid::beginCellCentroids ( ) const

Get an iterator over the cell centroids positioned at the first one.

◆ beginFaceCentroids()

CentroidIterator< 1 > Dune::CpGrid::beginFaceCentroids ( ) const

Get an iterator over the face centroids positioned at the first one.

◆ boundaryId()

int Dune::CpGrid::boundaryId ( int  face) const

◆ cellCenterDepth()

double Dune::CpGrid::cellCenterDepth ( int  cell_index) const

Get vertical position of cell center ("zcorn" average).

cell_index The index of the specific cell.

◆ cellCentroid()

const Vector & Dune::CpGrid::cellCentroid ( int  cell) const

Get the coordinates of the center of a cell.

Parameters
cellThe index identifying the face.

◆ cellCommunication()

const CommunicationType & Dune::CpGrid::cellCommunication ( ) const

Get the owner-overlap-copy communication for cells.

Suitable e.g. for parallel linear algebra used by CCFV

◆ cellFace()

int Dune::CpGrid::cellFace ( int  cell,
int  local_index 
) const

Get a specific face of a cell.

Parameters
cellThe index identifying the cell.
local_indexThe local index (in [0,numFaces(cell))) of the face in this cell.
Returns
The index identifying the face.

◆ cellFaceRow()

const cpgrid::OrientedEntityTable< 0, 1 >::row_type Dune::CpGrid::cellFaceRow ( int  cell) const

Get a list of indices identifying all faces of a cell.

Parameters
cellThe index identifying the cell.

◆ cellScatterGatherInterface()

const InterfaceMap & Dune::CpGrid::cellScatterGatherInterface ( ) const

Get an interface for gathering/scattering data attached to cells with communication.

Scattering means sending data from the indices of the global grid on process 0 to the distributed grid on all ranks independent of the grid. Gathering is the other way around. The interface can be used with VariableSizeCommunicator and a custom index based data handle to scatter (forward direction of the communicator) and gather data (backward direction of the communicator). Here is a small example that prints the received values when scattering:

struct Handle{
typedef int DataType;
const std::vector<int>& vals;
bool fixedsize() { return true; }
size_t size(std::size_t) { return 1; }
void gather(auto& B buf, size_t i)[ buf.write(vals[i]); }
void scatter(auto& B buf, size_t i, std::size_t) {
int val;
buf.read(val);
cout<<i<<": "<<val<<" "; }
};
Handle handle;
handle.vals.resize(grid.size(0), -1);
Dune::VariableSizeCommunicator<> comm(grid.comm(),
grid.cellScatterGatherInterface());
comm.forward(handle);
const CpGridTraits::Communication & comm() const
Get the collective communication object.
int size(int level, int codim) const
Number of grid entities per level and codim.

Referenced by scatterData().

◆ cellVolume()

double Dune::CpGrid::cellVolume ( int  cell) const

Get the volume of the cell.

Parameters
cellThe index identifying the cell.

◆ chooseData()

const std::vector< std::shared_ptr< Dune::cpgrid::CpGridData > > & Dune::CpGrid::chooseData ( ) const

◆ comm()

const CpGridTraits::Communication & Dune::CpGrid::comm ( ) const

Get the collective communication object.

◆ communicate() [1/2]

template<class DataHandle >
void Dune::CpGrid::communicate ( DataHandle &  data,
InterfaceType  iftype,
CommunicationDirection  dir 
) const

communicate objects for all codims on a given level.

The new communication interface.

Template Parameters
DataHandleThe type of the data handle describing the data.
Parameters
dataThe data handle describing the data. Has to adhere to the Dune::DataHandleIF interface.
iftypeThe interface to use for the communication.
dirThe direction of the communication along the interface (forward or backward).
leveldiscarded as CpGrid is not adaptive.

References Dune::cpgrid::CpGridData::communicate(), and data.

◆ communicate() [2/2]

template<class DataHandle >
void Dune::CpGrid::communicate ( DataHandle &  data,
InterfaceType  iftype,
CommunicationDirection  dir,
int   
) const
inline

communicate objects for all codims on a given level

The new communication interface.

Parameters
dataThe data handle describing the data. Has to adhere to the Dune::DataHandleIF interface.
iftypeThe interface to use for the communication.
dirThe direction of the communication along the interface (forward or backward).
leveldiscarded as CpGrid is not adaptive.

References communicate(), and data.

Referenced by communicate().

◆ createCartesian()

void Dune::CpGrid::createCartesian ( const std::array< int, 3 > &  dims,
const std::array< double, 3 > &  cellsize,
const std::array< int, 3 > &  shift = {0, 0, 0} 
)

Create a cartesian grid.

Parameters
dimsthe number of cells in each cartesian direction.
cellsizethe size of each cell in each dimension.
shiftThe origin of the grid, i.e. the corner of the cell with index (0,0,0) where the left, bottom, and top face of that cell intersect, is at the coordinate origin per default. This parameter shifts that corner to lie at (shift[0]*cellsize[0], ..., shift[2]*cellsize[2]).

◆ entity()

template<int codim>
cpgrid::Entity< codim > Dune::CpGrid::entity ( const cpgrid::Entity< codim > &  seed) const

given an EntitySeed (or EntityPointer) return an entity object

◆ faceArea()

double Dune::CpGrid::faceArea ( int  face) const

Get the area of a face.

Parameters
cellThe index identifying the face.

◆ faceAreaNormalEcl()

const Vector Dune::CpGrid::faceAreaNormalEcl ( int  face) const

◆ faceCell()

int Dune::CpGrid::faceCell ( int  face,
int  local_index 
) const

Get the index identifying a cell attached to a face.

Note that a face here is always oriented. If there are two neighboring cells then the orientation will be from local_index 0 to local_index 1

Parameters
faceThe index identifying the face.
local_indexThe local_index of the cell.
Returns
The index identifying a cell or -1 if there is no such cell due the face being part of the grid boundary or the cell being stored on another process.

Referenced by Dune::cpgrid::FaceCellsContainerProxy::operator()(), and Dune::cpgrid::FaceCellsProxy::operator[]().

◆ faceCenterEcl()

const Vector Dune::CpGrid::faceCenterEcl ( int  cell_index,
int  face,
const Dune::cpgrid::Intersection intersection 
) const

◆ faceCentroid()

const Vector & Dune::CpGrid::faceCentroid ( int  face) const

Get the coordinates of the center of a face.

Parameters
cellThe index identifying the face.

◆ faceNormal()

const Vector & Dune::CpGrid::faceNormal ( int  face) const

Get the unit normal of a face.

Parameters
cellThe index identifying the face.
See also
faceCell

◆ faceTag()

template<class Cell2FacesRowIterator >
int Dune::CpGrid::faceTag ( const Cell2FacesRowIterator &  cell_face) const

Get the cartesian tag associated with a face tag.

The tag tells us in which direction the face would point in the underlying cartesian grid.

Parameters
Aniterator that points to the face and was obtained by iterating over Opm::UgGridHelpers::cell2Faces(grid).

References I_FACE, J_FACE, K_FACE, NNC_FACE, numCells(), and numFaces().

◆ faceVertex()

int Dune::CpGrid::faceVertex ( int  face,
int  local_index 
) const

Get the index identifying a vertex of a face.

Parameters
cellThe index identifying the face.
local_indexThe local_index (in [0,numFaceVertices(vertex) - 1]]) of the vertex.

◆ gatherData()

template<class DataHandle >
void Dune::CpGrid::gatherData ( DataHandle &  handle) const

Moves data from the distributed view to the global (all data on process) view.

Template Parameters
DataHandleThe type of the data handle describing the data and responsible for gathering and scattering the data.
Parameters
handleThe data handle describing the data and responsible for gathering and scattering the data.

◆ geomTypes()

const std::vector< Dune::GeometryType > & Dune::CpGrid::geomTypes ( const int  ) const

◆ getCellIndexSet() [1/2]

ParallelIndexSet & Dune::CpGrid::getCellIndexSet ( )

◆ getCellIndexSet() [2/2]

const ParallelIndexSet & Dune::CpGrid::getCellIndexSet ( ) const

◆ getCellRemoteIndices() [1/2]

RemoteIndices & Dune::CpGrid::getCellRemoteIndices ( )

◆ getCellRemoteIndices() [2/2]

const RemoteIndices & Dune::CpGrid::getCellRemoteIndices ( ) const

◆ getEclCentroid() [1/2]

std::array< double, 3 > Dune::CpGrid::getEclCentroid ( const cpgrid::Entity< 0 > &  elem) const

◆ getEclCentroid() [2/2]

std::array< double, 3 > Dune::CpGrid::getEclCentroid ( const int &  idx) const

◆ getIJK()

void Dune::CpGrid::getIJK ( const int  c,
std::array< int, 3 > &  ijk 
) const

Extract Cartesian index triplet (i,j,k) of an active cell.

Parameters
[in]cActive cell index.
[out]ijkCartesian index triplet

Referenced by Dune::CartesianIndexMapper< CpGrid >::cartesianCoordinate().

◆ getLgrNameToLevel()

const std::map< std::string, int > & Dune::CpGrid::getLgrNameToLevel ( ) const

◆ getParentIntersectionFromLgrBoundaryFace()

Dune::cpgrid::Intersection Dune::CpGrid::getParentIntersectionFromLgrBoundaryFace ( const Dune::cpgrid::Intersection intersection) const

◆ ghostSize() [1/2]

unsigned int Dune::CpGrid::ghostSize ( int  ) const

Size of the ghost cell layer on the leaf level.

◆ ghostSize() [2/2]

unsigned int Dune::CpGrid::ghostSize ( int  ,
int   
) const

Size of the ghost cell layer on a given level.

◆ globalCell()

const std::vector< int > & Dune::CpGrid::globalCell ( ) const

Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells. Specifically,

const std::vector< int > & globalCell() const

is the linearized Cartesian index of grid cell

i

. This method should only be used by classes which really need it, such as those dealing with permeability fields from the input deck from whence the current CpGrid was constructed.

Referenced by Dune::CartesianIndexMapper< CpGrid >::cartesianIndex(), and Dune::CartesianIndexMapper< CpGrid >::compressedSize().

◆ globalIdSet()

const Traits::GlobalIdSet & Dune::CpGrid::globalIdSet ( ) const

Access to the GlobalIdSet.

◆ globalRefine()

void Dune::CpGrid::globalRefine ( int  )

global refinement

◆ lbegin() [1/2]

template<int codim>
Traits::template Codim< codim >::LevelIterator Dune::CpGrid::lbegin ( int  level) const

Iterator to first entity of given codim on level.

◆ lbegin() [2/2]

template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator Dune::CpGrid::lbegin ( int  level) const

Iterator to first entity of given codim on level and PartitionIteratorType.

◆ leafbegin() [1/2]

template<int codim>
Traits::template Codim< codim >::LeafIterator Dune::CpGrid::leafbegin ( ) const

Iterator to first leaf entity of given codim.

◆ leafbegin() [2/2]

template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator Dune::CpGrid::leafbegin ( ) const

Iterator to first leaf entity of given codim and PartitionIteratorType.

◆ leafend() [1/2]

template<int codim>
Traits::template Codim< codim >::LeafIterator Dune::CpGrid::leafend ( ) const

one past the end of the sequence of leaf entities

◆ leafend() [2/2]

template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator Dune::CpGrid::leafend ( ) const

one past the end of the sequence of leaf entities

◆ leafIndexSet()

const Traits::LeafIndexSet & Dune::CpGrid::leafIndexSet ( ) const

Access to the LeafIndexSet.

◆ lend() [1/2]

template<int codim>
Traits::template Codim< codim >::LevelIterator Dune::CpGrid::lend ( int  level) const

one past the end on this level

◆ lend() [2/2]

template<int codim, PartitionIteratorType PiType>
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator Dune::CpGrid::lend ( int  level) const

one past the end on this level

◆ levelIndexSet()

const Traits::LevelIndexSet & Dune::CpGrid::levelIndexSet ( int  level) const

Access to the LevelIndexSets.

◆ loadBalance() [1/9]

std::pair< bool, std::vector< std::pair< std::string, bool > > > Dune::CpGrid::loadBalance ( const std::vector< cpgrid::OpmWellType > *  wells,
const double *  transmissibilities = nullptr,
int  overlapLayers = 1,
bool  useZoltan = true 
)
inline

Distributes this grid over the available nodes in a distributed machine.

This will construct the corresponding graph to the grid and use the transmissibilities specified as weights associated with its edges. The graph will be passed to the load balancer.

Parameters
wellsThe wells of the eclipse If null wells will be neglected. If this is not null then complete well information of of the last scheduler step of the eclipse state will be used to make sure that all the possible completion cells of each well are stored on one process. This done by adding an edge with a very high edge weight for all possible pairs of cells in the completion set of a well.
transmissibilitiesThe transmissibilities used as the edge weights.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
Warning
May only be called once.
Returns
A pair consisting of a boolean indicating whether loadbalancing actually happened and a vector containing a pair of name and a boolean, indicating whether this well has perforated cells local to the process, for all wells (sorted by name)

References Dune::defaultTransEdgeWgt.

◆ loadBalance() [2/9]

bool Dune::CpGrid::loadBalance ( const std::vector< int > &  parts,
bool  ownersFirst = false,
bool  addCornerCells = false,
int  overlapLayers = 1 
)
inline

Distributes this grid over the available nodes in a distributed machine.

Parameters
partsThe partitioning information. For a cell with local index i the entry parts[i] is the partion number. Partition numbers need to start with zero and need to be consectutive also parts.size()==grid.leafGridView().size() and the ranks communicator need to be able to map all parts. Needs to valid at rank 0. Number of parts cannot exceed the number of ranks. Parts need to numbered consecutively starting from zero.
ownersFirstOrder owner cells before copy/overlap cells.
addCornerCellsAdd corner cells to the overlap layer.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
Warning
May only be called once.

References Dune::defaultTransEdgeWgt.

◆ loadBalance() [3/9]

template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > Dune::CpGrid::loadBalance ( DataHandle &  data,
const std::vector< cpgrid::OpmWellType > *  wells,
const double *  transmissibilities = nullptr,
int  overlapLayers = 1,
bool  useZoltan = true 
)
inline

Distributes this grid and data over the available nodes in a distributed machine.

Parameters
dataA data handle describing how to distribute attached data.
wellsThe wells of the eclipse Default: null If this is not null then complete well information of of the last scheduler step of the eclipse state will be used to make sure that all the possible completion cells of each well are stored on one process. This done by adding an edge with a very high edge weight for all possible pairs of cells in the completion set of a well.
transmissibilitiesThe transmissibilities used to calculate the edge weights.
overlapLayersThe number of layers of overlap cells to be added (default: 1)
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
Template Parameters
DataHandleThe type implementing DUNE's DataHandle interface.
Warning
May only be called once.
Returns
A pair consisting of a boolean indicating whether loadbalancing actually happened and a vector containing a pair of name and a boolean, indicating whether this well has perforated cells local to the process, for all wells (sorted by name)

References data, loadBalance(), and scatterData().

◆ loadBalance() [4/9]

template<class DataHandle >
bool Dune::CpGrid::loadBalance ( DataHandle &  data,
const std::vector< int > &  parts,
bool  ownersFirst = false,
bool  addCornerCells = false,
int  overlapLayers = 1 
)
inline

Distributes this grid and data over the available nodes in a distributed machine.

Parameters
dataA data handle describing how to distribute attached data.
partsThe partitioning information. For a cell with local index i the entry parts[i] is the partion number. Partition numbers need to start with zero and need to be consectutive also parts.size()==grid.leafGridView().size() and the ranks communicator need to be able to map all parts. Needs to valid at rank 0. Number of parts cannot exceed the number of ranks. Parts need to numbered consecutively starting from zero.
ownersFirstOrder owner cells before copy/overlap cells.
addCornerCellsAdd corner cells to the overlap layer.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
Warning
May only be called once.

References data, loadBalance(), and scatterData().

◆ loadBalance() [5/9]

template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > Dune::CpGrid::loadBalance ( DataHandle &  data,
const std::vector< int > &  parts,
const std::vector< cpgrid::OpmWellType > *  wells,
bool  ownersFirst = false,
bool  addCornerCells = false,
int  overlapLayers = 1 
)
inline

Distributes this grid over the available nodes in a distributed machine.

Parameters
dataA data handle describing how to distribute attached data.
partsThe partitioning information. For a cell with local index i the entry parts[i] is the partion number. Partition numbers need to start with zero and need to be consectutive also parts.size()==grid.leafGridView().size() and the ranks communicator need to be able to map all parts. Needs to valid at rank 0. Number of parts cannot exceed the number of ranks. Parts need to numbered consecutively starting from zero.
ownersFirstOrder owner cells before copy/overlap cells.
addCornerCellsAdd corner cells to the overlap layer.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
Template Parameters
DataHandleThe type implementing DUNE's DataHandle interface.
Warning
May only be called once.
Returns
A pair consisting of a boolean indicating whether loadbalancing actually happened and a vector containing a pair of name and a boolean, indicating whether this well has perforated cells local to the process, for all wells (sorted by name)

References data, Dune::defaultTransEdgeWgt, and scatterData().

◆ loadBalance() [6/9]

template<class DataHandle >
bool Dune::CpGrid::loadBalance ( DataHandle &  data,
decltype(data.fixedSize(0, 0))  overlapLayers = 1,
bool  useZoltan = true 
)
inline

Distributes this grid and data over the available nodes in a distributed machine.

Parameters
dataA data handle describing how to distribute attached data.
overlapLayersThe number of layers of overlap cells to be added (default: 1)
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
Template Parameters
DataHandleThe type implementing DUNE's DataHandle interface.
Warning
May only be called once.

References data, loadBalance(), and scatterData().

◆ loadBalance() [7/9]

template<class DataHandle >
std::pair< bool, std::vector< std::pair< std::string, bool > > > Dune::CpGrid::loadBalance ( DataHandle &  data,
EdgeWeightMethod  method,
const std::vector< cpgrid::OpmWellType > *  wells,
bool  serialPartitioning,
const double *  transmissibilities = nullptr,
bool  ownersFirst = false,
bool  addCornerCells = false,
int  overlapLayers = 1,
bool  useZoltan = true,
double  zoltanImbalanceTol = 1.1,
bool  allowDistributedWells = false 
)
inline

Distributes this grid over the available nodes in a distributed machine.

This will construct the corresponding graph to the grid and use the transmissibilities specified to calculate the weights associated with its edges. The graph will be passed to the load balancer.

Parameters
dataA data handle describing how to distribute attached data.
methodThe edge-weighting method to be used on the Zoltan partitioner.
wellsThe information about all possible wells. If null then the wells will be neglected. Otherwise the wells will be used to make sure that all the possible completion cells of each well are stored on one process. This is done by adding an edge with a very high edge weight for all possible pairs of cells in the completion set of a well.
serialPartitioningIf true, the partitioning will be done on a single process.
transmissibilitiesThe transmissibilities used to calculate the edge weights.
ownersFirstOrder owner cells before copy/overlap cells.
addCornerCellsAdd corner cells to the overlap layer.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
zoltanImbalanceTolSet the imbalance tolerance used by Zoltan
allowDistributedWellsAllow the perforation of a well to be distributed to the interior region of multiple processes.
Template Parameters
DataHandleThe type implementing DUNE's DataHandle interface.
Warning
May only be called once.
Returns
A pair consisting of a boolean indicating whether loadbalancing actually happened and a vector containing a pair of name and a boolean, indicating whether this well has perforated cells local to the process, for all wells (sorted by name)

References data, and scatterData().

◆ loadBalance() [8/9]

std::pair< bool, std::vector< std::pair< std::string, bool > > > Dune::CpGrid::loadBalance ( EdgeWeightMethod  method,
const std::vector< cpgrid::OpmWellType > *  wells,
const double *  transmissibilities = nullptr,
bool  ownersFirst = false,
bool  addCornerCells = false,
int  overlapLayers = 1,
bool  useZoltan = true 
)
inline

Distributes this grid over the available nodes in a distributed machine.

This will construct the corresponding graph to the grid and use the transmissibilities specified to calculate the weights associated with its edges. The graph will be passed to the load balancer.

Parameters
methodThe edge-weighting method to be used on the Zoltan partitioner.
wellsThe wells of the eclipse If null wells will be neglected. If this is not null then complete well information of of the last scheduler step of the eclipse state will be used to make sure that all the possible completion cells of each well are stored on one process. This done by adding an edge with a very high edge weight for all possible pairs of cells in the completion set of a well.
transmissibilitiesThe transmissibilities used to calculate the edge weights.
ownersFirstOrder owner cells before copy/overlap cells.
addCornerCellsAdd corner cells to the overlap layer.
overlapLayersThe number of layers of cells of the overlap region (default: 1).
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
Warning
May only be called once.
Returns
A pair consisting of a boolean indicating whether loadbalancing actually happened and a vector containing a pair of name and a boolean, indicating whether this well has perforated cells local to the process, for all wells (sorted by name)

◆ loadBalance() [9/9]

bool Dune::CpGrid::loadBalance ( int  overlapLayers = 1,
bool  useZoltan = true 
)
inline

Distributes this grid over the available nodes in a distributed machine.

Parameters
overlapLayersThe number of layers of cells of the overlap region (default: 1).
useZoltanWhether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid.
Warning
May only be called once.

References Dune::defaultTransEdgeWgt.

Referenced by loadBalance().

◆ localIdSet()

const Traits::LocalIdSet & Dune::CpGrid::localIdSet ( ) const

Access to the LocalIdSet.

◆ logicalCartesianSize()

const std::array< int, 3 > & Dune::CpGrid::logicalCartesianSize ( ) const

The logical cartesian size of the global grid. This function is not part of the Dune grid interface, and should be used with caution.

Referenced by Dune::CartesianIndexMapper< CpGrid >::cartesianDimensions().

◆ maxLevel()

int Dune::CpGrid::maxLevel ( ) const

Return maximum level defined in this grid. Levels are 0 and 1, maxlevel = 1 (not counting leafview), 0 = the coarsest level.

Referenced by Dune::CartesianIndexMapper< CpGrid >::cartesianCoordinateLevel().

◆ name()

std::string Dune::CpGrid::name ( ) const

Get the grid name.

It's the same as the class name. What did you expect, something funny?

◆ numBoundarySegments()

unsigned int Dune::CpGrid::numBoundarySegments ( ) const

returns the number of boundary segments within the macro grid

◆ numCellFaces() [1/2]

int Dune::CpGrid::numCellFaces ( ) const

Get the sum of all faces attached to all cells.

Each face identified by a unique index is counted as often as there are neigboring cells attached to it. $ numCellFaces()=\sum_{c} numCellFaces(c) $

See also
numCellFaces(int)const

◆ numCellFaces() [2/2]

int Dune::CpGrid::numCellFaces ( int  cell) const

Get the number of faces of a cell.

Due to faults, and collapsing vertices (along pillars) this number is quite arbitrary. Its lower bound is 4, but there is no upper bound. \parame cell the index identifying the cell.

◆ numCells()

int Dune::CpGrid::numCells ( ) const

Get the number of cells.

Referenced by faceTag().

◆ numFaces()

int Dune::CpGrid::numFaces ( ) const

Get the number of faces.

Referenced by faceTag().

◆ numFaceVertices()

int Dune::CpGrid::numFaceVertices ( int  face) const

◆ numVertices()

int Dune::CpGrid::numVertices ( ) const

Get The number of vertices.

◆ overlapSize() [1/2]

unsigned int Dune::CpGrid::overlapSize ( int  ) const

Size of the overlap on the leaf level.

◆ overlapSize() [2/2]

unsigned int Dune::CpGrid::overlapSize ( int  ,
int   
) const

Size of the overlap on a given level.

◆ pointScatterGatherInterface()

const InterfaceMap & Dune::CpGrid::pointScatterGatherInterface ( ) const

Get an interface for gathering/scattering data attached to points with communication.

See also
cellScatterGatherInterface

Referenced by scatterData().

◆ processEclipseFormat()

void Dune::CpGrid::processEclipseFormat ( const grdecl input_data,
bool  remove_ij_boundary,
bool  turn_normals = false 
)

Read the Eclipse grid format ('grdecl').

Parameters
input_datathe data in grdecl format, declared in preprocess.h.
remove_ij_boundaryif true, will remove (i, j) boundaries. Used internally.

◆ readSintefLegacyFormat()

void Dune::CpGrid::readSintefLegacyFormat ( const std::string &  grid_prefix)

Read the Sintef legacy grid format ('topogeom').

Parameters
grid_prefixthe grid name, such that topology is found in <grid_prefix>-topo.dat etc.

◆ scatterData()

template<class DataHandle >
void Dune::CpGrid::scatterData ( DataHandle &  handle) const

Moves data from the global (all data on process) view to the distributed view.

This method does not do communication but assumes that the global grid is present on every process and simply copies data to the distributed view.

Template Parameters
DataHandleThe type of the data handle describing the data and responsible for gathering and scattering the data.
Parameters
handleThe data handle describing the data and responsible for gathering and scattering the data.

References cellScatterGatherInterface(), and pointScatterGatherInterface().

Referenced by loadBalance().

◆ setUniqueBoundaryIds()

void Dune::CpGrid::setUniqueBoundaryIds ( bool  uids)

Set whether we want to have unique boundary ids.

Parameters
uidsif true, each boundary intersection will have a unique boundary id.

◆ setZoltanParams()

void Dune::CpGrid::setZoltanParams ( const std::map< std::string, std::string > &  params)

◆ size() [1/4]

int Dune::CpGrid::size ( GeometryType  type) const

number of leaf entities per geometry type in this process

◆ size() [2/4]

int Dune::CpGrid::size ( int  codim) const

number of leaf entities per codim in this process

◆ size() [3/4]

int Dune::CpGrid::size ( int  level,
GeometryType  type 
) const

number of entities per level and geometry type in this process

◆ size() [4/4]

int Dune::CpGrid::size ( int  level,
int  codim 
) const

Number of grid entities per level and codim.

◆ sortedNumAquiferCells()

const std::vector< int > & Dune::CpGrid::sortedNumAquiferCells ( ) const

Get sorted active cell indices of numerical aquifer.

◆ switchToDistributedView()

void Dune::CpGrid::switchToDistributedView ( )

Switch to the distributed view.

◆ switchToGlobalView()

void Dune::CpGrid::switchToGlobalView ( )

Switch to the global view.

◆ uniqueBoundaryIds()

bool Dune::CpGrid::uniqueBoundaryIds ( ) const

Is the grid currently using unique boundary ids?

Returns
true if each boundary intersection has a unique id false if we use the (default) 1-6 ids for i- i+ j- j+ k- k+ boundaries.

◆ vertexPosition()

const Vector & Dune::CpGrid::vertexPosition ( int  vertex) const

Get the Position of a vertex.

Parameters
cellThe index identifying the cell.
Returns
The coordinates of the vertex.

◆ writeSintefLegacyFormat()

void Dune::CpGrid::writeSintefLegacyFormat ( const std::string &  grid_prefix) const

Write the Sintef legacy grid format ('topogeom').

Parameters
grid_prefixthe grid name, such that topology will be found in <grid_prefix>-topo.dat etc.

◆ zcornData()

const std::vector< double > & Dune::CpGrid::zcornData ( ) const

◆ zoltanPartitionWithoutScatter()

std::vector< int > Dune::CpGrid::zoltanPartitionWithoutScatter ( const std::vector< cpgrid::OpmWellType > *  wells,
const double *  transmissibilities,
const int  numParts,
const double  zoltanImbalanceTol 
) const

Partitions the grid using Zoltan without decomposing and distributing it among processes.

Parameters
wellsThe wells of the eclipse.
transmissibilitiesThe transmissibilities used to calculate the edge weights.
numPartsNumber of parts in the partition.
Returns
An array with the domain index for each cell.

Friends And Related Function Documentation

◆ cpgrid::CpGridData

friend class cpgrid::CpGridData
friend

◆ cpgrid::Entity< 0 >

friend class cpgrid::Entity< 0 >
friend

◆ cpgrid::Entity< 1 >

friend class cpgrid::Entity< 1 >
friend

◆ cpgrid::Entity< 2 >

friend class cpgrid::Entity< 2 >
friend

◆ cpgrid::Entity< 3 >

friend class cpgrid::Entity< 3 >
friend

◆ createEntity

template<int dim>
cpgrid::Entity< dim > createEntity ( const CpGrid ,
int  ,
bool   
)
friend

◆ Opm::LookUpCartesianData

template<typename Grid , typename GridView >
friend class Opm::LookUpCartesianData
friend

◆ Opm::LookUpData

template<typename Grid , typename GridView >
friend class Opm::LookUpData
friend

The documentation for this class was generated from the following file: