[ provides Dune::Grid ]
More...
#include <CpGrid.hpp>
|
| CpGrid () |
| Default constructor. More...
|
|
| CpGrid (MPIHelper::MPICommunicator comm) |
|
|
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) |
|
|
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) |
|
|
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 > ¶ms) |
|
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::Communication & | comm () const |
| Get the collective communication object. More...
|
|
◆ 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
The type of the parallel index set.
◆ RemoteIndices
The type of the remote indices information.
◆ Vector
◆ CpGrid() [1/2]
◆ CpGrid() [2/2]
Dune::CpGrid::CpGrid |
( |
MPIHelper::MPICommunicator |
comm | ) |
|
◆ 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_vec | Vector of Number of (refined) cells in each direction that each parent cell should be refined to. |
[in] | startIJK_vec | Vector of Cartesian triplet indices where each patch starts. |
[in] | endIJK_vec | Vector 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_vec | Names (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_dim | Number of (refined) cells in each direction that each parent cell should be refined to. |
[in] | startIJK | Cartesian triplet index where the patch starts. |
[in] | endIJK | Cartesian triplet index where the patch ends. Last cell part of the lgr will be {endijk[0]-1, ... endIJK[2]-1}. |
[in] | lgr_name | Name (std::string) for the lgr/level1 |
◆ beginCellCentroids()
Get an iterator over the cell centroids positioned at the first one.
◆ beginFaceCentroids()
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
-
cell | The index identifying the face. |
◆ cellCommunication()
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
-
cell | The index identifying the cell. |
local_index | The local index (in [0,numFaces(cell))) of the face in this cell. |
- Returns
- The index identifying the face.
◆ cellFaceRow()
Get a list of indices identifying all faces of a cell.
- Parameters
-
cell | The 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());
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
-
cell | The index identifying the cell. |
◆ chooseData()
◆ comm()
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
-
DataHandle | The type of the data handle describing the data. |
- Parameters
-
data | The data handle describing the data. Has to adhere to the Dune::DataHandleIF interface. |
iftype | The interface to use for the communication. |
dir | The direction of the communication along the interface (forward or backward). |
level | discarded 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
-
data | The data handle describing the data. Has to adhere to the Dune::DataHandleIF interface. |
iftype | The interface to use for the communication. |
dir | The direction of the communication along the interface (forward or backward). |
level | discarded 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
-
dims | the number of cells in each cartesian direction. |
cellsize | the size of each cell in each dimension. |
shift | The 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()
given an EntitySeed (or EntityPointer) return an entity object
◆ faceArea()
double Dune::CpGrid::faceArea |
( |
int |
face | ) |
const |
Get the area of a face.
- Parameters
-
cell | The 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
-
face | The index identifying the face. |
local_index | The 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()
◆ faceCentroid()
const Vector & Dune::CpGrid::faceCentroid |
( |
int |
face | ) |
const |
Get the coordinates of the center of a face.
- Parameters
-
cell | The index identifying the face. |
◆ faceNormal()
const Vector & Dune::CpGrid::faceNormal |
( |
int |
face | ) |
const |
Get the unit normal of a face.
- Parameters
-
cell | The 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
-
An | iterator 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
-
cell | The index identifying the face. |
local_index | The 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
-
DataHandle | The type of the data handle describing the data and responsible for gathering and scattering the data. |
- Parameters
-
handle | The 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]
◆ getCellIndexSet() [2/2]
◆ getCellRemoteIndices() [1/2]
◆ 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 |
◆ getLgrNameToLevel()
const std::map< std::string, int > & Dune::CpGrid::getLgrNameToLevel |
( |
| ) |
const |
◆ getParentIntersectionFromLgrBoundaryFace()
◆ 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 |
◆ globalIdSet()
const Traits::GlobalIdSet & Dune::CpGrid::globalIdSet |
( |
| ) |
const |
Access to the GlobalIdSet.
◆ globalRefine()
void Dune::CpGrid::globalRefine |
( |
int |
| ) |
|
◆ 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
-
wells | The 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. |
transmissibilities | The transmissibilities used as the edge weights. |
overlapLayers | The number of layers of cells of the overlap region (default: 1). |
useZoltan | Whether 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
-
parts | The 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. |
ownersFirst | Order owner cells before copy/overlap cells. |
addCornerCells | Add corner cells to the overlap layer. |
overlapLayers | The 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
-
data | A data handle describing how to distribute attached data. |
wells | The 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. |
transmissibilities | The transmissibilities used to calculate the edge weights. |
overlapLayers | The number of layers of overlap cells to be added (default: 1) |
useZoltan | Whether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid. |
- Template Parameters
-
DataHandle | The 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
-
data | A data handle describing how to distribute attached data. |
parts | The 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. |
ownersFirst | Order owner cells before copy/overlap cells. |
addCornerCells | Add corner cells to the overlap layer. |
overlapLayers | The 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
-
data | A data handle describing how to distribute attached data. |
parts | The 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. |
ownersFirst | Order owner cells before copy/overlap cells. |
addCornerCells | Add corner cells to the overlap layer. |
overlapLayers | The number of layers of cells of the overlap region (default: 1). |
- Template Parameters
-
DataHandle | The 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
-
data | A data handle describing how to distribute attached data. |
overlapLayers | The number of layers of overlap cells to be added (default: 1) |
useZoltan | Whether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid. |
- Template Parameters
-
DataHandle | The 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
-
data | A data handle describing how to distribute attached data. |
method | The edge-weighting method to be used on the Zoltan partitioner. |
wells | The 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. |
serialPartitioning | If true, the partitioning will be done on a single process. |
transmissibilities | The transmissibilities used to calculate the edge weights. |
ownersFirst | Order owner cells before copy/overlap cells. |
addCornerCells | Add corner cells to the overlap layer. |
overlapLayers | The number of layers of cells of the overlap region (default: 1). |
useZoltan | Whether to use Zoltan for partitioning or our simple approach based on rectangular partitioning the underlying cartesian grid. |
zoltanImbalanceTol | Set the imbalance tolerance used by Zoltan |
allowDistributedWells | Allow the perforation of a well to be distributed to the interior region of multiple processes. |
- Template Parameters
-
DataHandle | The 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
-
method | The edge-weighting method to be used on the Zoltan partitioner. |
wells | The 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. |
transmissibilities | The transmissibilities used to calculate the edge weights. |
ownersFirst | Order owner cells before copy/overlap cells. |
addCornerCells | Add corner cells to the overlap layer. |
overlapLayers | The number of layers of cells of the overlap region (default: 1). |
useZoltan | Whether 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
-
overlapLayers | The number of layers of cells of the overlap region (default: 1). |
useZoltan | Whether 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 |
◆ maxLevel()
int Dune::CpGrid::maxLevel |
( |
| ) |
const |
◆ 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. - 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 |
◆ 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_data | the data in grdecl format, declared in preprocess.h. |
remove_ij_boundary | if 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_prefix | the 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
-
DataHandle | The type of the data handle describing the data and responsible for gathering and scattering the data. |
- Parameters
-
handle | The 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
-
uids | if 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
-
cell | The 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_prefix | the 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
-
wells | The wells of the eclipse. |
transmissibilities | The transmissibilities used to calculate the edge weights. |
numParts | Number of parts in the partition. |
- Returns
- An array with the domain index for each cell.
◆ cpgrid::CpGridData
◆ cpgrid::Entity< 0 >
◆ cpgrid::Entity< 1 >
◆ cpgrid::Entity< 2 >
◆ cpgrid::Entity< 3 >
◆ createEntity
◆ Opm::LookUpCartesianData
template<typename Grid , typename GridView >
◆ Opm::LookUpData
template<typename Grid , typename GridView >
The documentation for this class was generated from the following file:
|