38 #ifndef OPM_CPGRID_HEADER
39 #define OPM_CPGRID_HEADER
44 #include <unordered_set>
45 #include <opm/grid/utility/ErrorMacros.hpp>
48 #include <opm/grid/utility/platform_dependent/disable_warnings.h>
50 #include <dune/common/version.hh>
53 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
54 #include <dune/common/parallel/variablesizecommunicator.hh>
56 #include <opm/grid/utility/VariableSizeCommunicator.hpp>
60 #include <dune/grid/common/capabilities.hh>
61 #include <dune/grid/common/grid.hh>
62 #include <dune/grid/common/gridenums.hh>
64 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
66 #include "cpgrid/Intersection.hpp"
67 #include "cpgrid/Entity.hpp"
68 #include "cpgrid/Geometry.hpp"
69 #include "cpgrid/Iterators.hpp"
70 #include "cpgrid/Indexsets.hpp"
71 #include "cpgrid/DefaultGeometryPolicy.hpp"
72 #include "common/GridEnums.hpp"
73 #include "common/Volumes.hpp"
76 #include <opm/grid/utility/OpmParserIncludes.hpp>
145 template <PartitionIteratorType pitype>
157 template <PartitionIteratorType pitype>
163 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> >
LeafGridView;
170 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> >
LeafGridView;
184 typedef Dune::CollectiveCommunication<MPICommunicator> CollectiveCommunication;
206 :
public GridDefaultImplementation<3, 3, double, CpGridFamily >
263 Opm::EclipseState* ecl_state,
264 bool periodic_extension,
bool turn_normals,
bool clip_z,
287 Opm::EclipseState* ecl_state,
288 bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false);
309 const std::array<double, 3>& cellsize);
316 return current_view_data_->logical_cartesian_size_;
328 return current_view_data_->global_cell_;
338 void getIJK(
const int c, std::array<int,3>& ijk)
const
340 current_view_data_->
getIJK(c, ijk);
383 typename Traits::template Codim<codim>::LevelIterator
lbegin (
int level)
const
386 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
393 typename Traits::template Codim<codim>::LevelIterator
lend (
int level)
const
396 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
403 template<
int codim, PartitionIteratorType PiType>
404 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator
lbegin (
int level)
const
407 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
413 template<
int codim, PartitionIteratorType PiType>
414 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator
lend (
int level)
const
417 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
424 typename Traits::template Codim<codim>::LeafIterator
leafbegin()
const
432 typename Traits::template Codim<codim>::LeafIterator
leafend()
const
439 template<
int codim, PartitionIteratorType PiType>
440 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator
leafbegin()
const
447 template<
int codim, PartitionIteratorType PiType>
448 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator
leafend()
const
455 int size (
int level,
int codim)
const
458 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
466 return current_view_data_->
size(codim);
471 int size (
int level, GeometryType type)
const
474 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
480 int size (GeometryType type)
const
482 return current_view_data_->
size(type);
489 return global_id_set_;
496 return global_id_set_;
504 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
505 return *current_view_data_->index_set_;
512 return *current_view_data_->index_set_;
519 std::cout <<
"Warning: Global refinement not implemented, yet." << std::endl;
522 const std::vector< Dune :: GeometryType >& geomTypes(
const int codim )
const
607 return current_view_data_->unique_boundary_ids_.size();
611 unsigned int numBndSegs = 0;
613 for (
int i = 0; i < num_faces; ++i) {
615 if (current_view_data_->face_to_cell_[face].
size() == 1) {
634 return get<0>(scatterGrid(
defaultTransEdgeWgt,
false,
nullptr,
false,
nullptr,
true, overlapLayers, useZoltan ));
658 std::pair<bool, std::vector<std::pair<std::string,bool> > >
660 const double* transmissibilities =
nullptr,
661 int overlapLayers=1,
bool useZoltan=
true)
663 return scatterGrid(
defaultTransEdgeWgt,
false, wells,
false, transmissibilities,
false, overlapLayers, useZoltan);
691 std::pair<bool, std::vector<std::pair<std::string,bool> > >
693 const double* transmissibilities =
nullptr,
bool ownersFirst=
false,
694 bool addCornerCells=
false,
int overlapLayers=1,
695 bool useZoltan =
true)
697 return scatterGrid(method, ownersFirst, wells,
false, transmissibilities, addCornerCells, overlapLayers, useZoltan);
719 template<
class DataHandle>
720 std::pair<bool, std::vector<std::pair<std::string,bool> > >
722 const std::vector<cpgrid::OpmWellType> * wells,
723 const double* transmissibilities =
nullptr,
724 int overlapLayers=1,
bool useZoltan =
true)
726 auto ret =
loadBalance(wells, transmissibilities, overlapLayers, useZoltan);
763 template<
class DataHandle>
764 std::pair<bool, std::vector<std::pair<std::string,bool> > >
766 const std::vector<cpgrid::OpmWellType> * wells,
767 bool serialPartitioning,
768 const double* transmissibilities =
nullptr,
bool ownersFirst=
false,
769 bool addCornerCells=
false,
int overlapLayers=1,
bool useZoltan =
true,
770 double zoltanImbalanceTol = 1.1,
771 bool allowDistributedWells =
false)
773 auto ret = scatterGrid(method, ownersFirst, wells, serialPartitioning, transmissibilities,
774 addCornerCells, overlapLayers, useZoltan, zoltanImbalanceTol, allowDistributedWells);
799 template<
class DataHandle>
800 std::pair<bool, std::vector<std::pair<std::string,bool> > >
802 const std::vector<cpgrid::OpmWellType> * wells,
803 bool ownersFirst=
false,
804 bool addCornerCells=
false,
int overlapLayers=1)
810 addCornerCells, overlapLayers,
false,
828 template<
class DataHandle>
830 #
if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7)
831 decltype(data.fixedSize(0,0)) overlapLayers=1,
bool useZoltan =
true)
833 decltype(data.fixedsize(0,0)) overlapLayers=1,
bool useZoltan = true)
857 bool loadBalance(
const std::vector<int>& parts,
bool ownersFirst=
false,
858 bool addCornerCells=
false,
int overlapLayers=1)
864 addCornerCells, overlapLayers,
false,
881 template<
class DataHandle>
882 bool loadBalance(DataHandle& data,
const std::vector<int>& parts,
bool ownersFirst=
false,
883 bool addCornerCells=
false,
int overlapLayers=1)
885 bool ret =
loadBalance(parts, ownersFirst, addCornerCells, overlapLayers);
900 template<
class DataHandle>
901 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir,
int )
const
913 template<
class DataHandle>
914 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir)
const
916 current_view_data_->
communicate(data, iftype, dir);
920 const CollectiveCommunication&
comm ()
const
922 return current_view_data_->ccobj_;
936 typedef Dune::FieldVector<double, 3> Vector;
939 const std::vector<double>& zcornData()
const {
940 return data_->zcornData();
948 return current_view_data_->cell_to_face_.
size();
953 return current_view_data_->face_to_cell_.
size();
958 return current_view_data_->geomVector<3>().
size();
976 return current_view_data_->cell_to_face_[
cpgrid::EntityRep<0>(cell,
true)][local_index].index();
1002 bool a = (local_index == 0);
1003 bool b = r[0].orientation();
1004 bool use_first = a ? b : !b;
1006 int r_size = r.size();
1009 if(r[0].index()==std::numeric_limits<int>::max()){
1014 if(r.size()>1 && r[1].index()==std::numeric_limits<int>::max())
1020 return use_first ? r[0].index() : r[1].index();
1022 return use_first ? r[index].index() : -1;
1033 return current_view_data_->cell_to_face_.
dataSize();
1035 int numFaceVertices(
int face)
const
1037 return current_view_data_->face_to_point_[face].
size();
1045 return current_view_data_->face_to_point_[face][local_index];
1054 const int nv = current_view_data_->cell_to_point_[cell_index].size();
1056 for (
int i=0; i<nv; ++i) {
1057 zz +=
vertexPosition(current_view_data_->cell_to_point_[cell_index][i])[nd-1];
1062 const Vector faceCenterEcl(
int cell_index,
int face)
const
1077 static const int faceVxMap[ 6 ][ 4 ] = { {0, 2, 4, 6},
1086 assert (current_view_data_->cell_to_point_[cell_index].size() == 8);
1088 for(
int i=0; i<4; ++i )
1090 center +=
vertexPosition(current_view_data_->cell_to_point_[cell_index][ faceVxMap[ face ][ i ] ]);
1093 for (
int i=0; i<3; ++i) {
1100 const Vector faceAreaNormalEcl(
int face)
const
1103 const int nd = Vector::dimension;
1104 const int nv = numFaceVertices(face);
1118 Vector areaNormal =
cross(a,b);
1119 for (
int i=0; i<nd; ++i) {
1129 Vector areaNormal =
cross(a,b);
1137 int k = (nv % 2) ? 0 : nv - 1;
1139 Vector areaNormal(0.0);
1141 for (
int i = 1; i < h; ++i)
1143 Vector a =
vertexPosition(current_view_data_->face_to_point_[face][2*i]) -
vertexPosition(current_view_data_->face_to_point_[face][0]);
1144 Vector b =
vertexPosition(current_view_data_->face_to_point_[face][2*i+1]) -
vertexPosition(current_view_data_->face_to_point_[face][2*i-1]);
1145 areaNormal +=
cross(a,b);
1149 Vector a =
vertexPosition(current_view_data_->face_to_point_[face][2*h]) -
vertexPosition(current_view_data_->face_to_point_[face][0]);
1150 Vector b =
vertexPosition(current_view_data_->face_to_point_[face][k]) -
vertexPosition(current_view_data_->face_to_point_[face][2*h-1]);
1151 areaNormal +=
cross(a,b);
1186 return current_view_data_->face_normals_.get(face);
1205 :
public RandomAccessIteratorFacade<CentroidIterator<codim>,
1206 FieldVector<double, 3>,
1207 const FieldVector<double, 3>&, int>
1211 typedef typename std::vector<
cpgrid::Geometry<3-codim, 3> >::const_iterator
1219 const FieldVector<double, 3>& dereference()
const
1221 return iter_->center();
1227 const FieldVector<double, 3>& elementAt(
int n)
1229 return iter_[n]->center();
1265 int boundaryId(
int face)
const
1274 if (current_view_data_->face_to_cell_[f].
size() == 1) {
1277 ret = current_view_data_->unique_boundary_ids_[f];
1280 const bool normal_is_in =
1281 !(current_view_data_->face_to_cell_[f][0].orientation());
1282 enum face_tag tag = current_view_data_->face_tag_[f];
1286 ret = normal_is_in ? 1 : 2;
1290 ret = normal_is_in ? 3 : 4;
1295 ret = normal_is_in ? 5 : 6;
1300 OPM_THROW(std::logic_error,
"NNC face at boundary. This should never happen!");
1313 template<
class Cell2FacesRowIterator>
1315 faceTag(
const Cell2FacesRowIterator& cell_face)
const
1326 const int cell = cell_face.getCellIndex();
1327 const int face = *cell_face;
1328 assert (0 <= cell); assert (cell <
numCells());
1329 assert (0 <= face); assert (face <
numFaces());
1334 const F2C& f2c = current_view_data_->face_to_cell_[f];
1335 const face_tag tag = current_view_data_->face_tag_[f];
1337 assert ((f2c.size() == 1) || (f2c.size() == 2));
1339 int inside_cell = 0;
1341 if ( f2c.size() == 2 )
1343 if ( f2c[1].index() == cell )
1348 const bool normal_is_in = ! f2c[inside_cell].orientation();
1353 return normal_is_in ? 0 : 1;
1356 return normal_is_in ? 2 : 3;
1360 return normal_is_in ? 4 : 5;
1365 OPM_THROW(std::logic_error,
"Unhandled face tag. This should never happen!");
1380 template<
class DataHandle>
1393 if(!distributed_data_)
1394 OPM_THROW(std::runtime_error,
"Moving Data only allowed with a load balanced grid!");
1409 template<
class DataHandle>
1413 if(!distributed_data_)
1414 OPM_THROW(std::runtime_error,
"Moving Data only allowed with a load balance grid!");
1415 distributed_data_->gatherData(handle, data_.get(), distributed_data_.get());
1422 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 7)
1424 using InterfaceMap = VariableSizeCommunicator<>::InterfaceMap;
1427 using InterfaceMap = Opm::VariableSizeCommunicator<>::InterfaceMap;
1467 return *cell_scatter_gather_interfaces_;
1474 return *point_scatter_gather_interfaces_;
1480 current_view_data_=data_.get();
1486 if (! distributed_data_)
1487 OPM_THROW(std::logic_error,
"No distributed view available in grid");
1488 current_view_data_=distributed_data_.get();
1494 typedef cpgrid::CpGridData::ParallelIndexSet ParallelIndexSet;
1496 typedef cpgrid::CpGridData::RemoteIndices RemoteIndices;
1499 using CommunicationType = cpgrid::CpGridData::CommunicationType;
1504 const CommunicationType& cellCommunication()
const
1506 return current_view_data_->cellCommunication();
1509 ParallelIndexSet& getCellIndexSet()
1511 return current_view_data_->cellIndexSet();
1514 RemoteIndices& getCellRemoteIndices()
1516 return current_view_data_->cellRemoteIndices();
1519 const ParallelIndexSet& getCellIndexSet()
const
1521 return current_view_data_->cellIndexSet();
1524 const RemoteIndices& getCellRemoteIndices()
const
1526 return current_view_data_->cellRemoteIndices();
1557 std::pair<bool, std::vector<std::pair<std::string,bool> > >
1560 const std::vector<cpgrid::OpmWellType> * wells,
1561 bool serialPartitioning,
1562 const double* transmissibilities,
1563 bool addCornerCells,
1565 bool useZoltan =
true,
1566 double zoltanImbalanceTol = 1.1,
1567 bool allowDistributedWells =
true,
1568 const std::vector<int>& input_cell_part = {});
1574 std::shared_ptr<cpgrid::CpGridData> data_;
1576 cpgrid::CpGridData* current_view_data_;
1578 std::shared_ptr<cpgrid::CpGridData> distributed_data_;
1584 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1590 std::shared_ptr<InterfaceMap> point_scatter_gather_interfaces_;
1594 cpgrid::GlobalIdSet global_id_set_;
1599 namespace Capabilities
1605 static const bool v =
true;
1612 static const bool v =
true;
1618 static const bool v =
true;
1624 static const bool v =
true;
1631 static const bool v =
false;
1645 #include <opm/grid/cpgrid/PersistentContainer.hpp>
1646 #include <opm/grid/cpgrid/CartesianIndexMapper.hpp>
An iterator over the centroids of the geometry of the entities.
Definition: CpGrid.hpp:1208
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition: CpGrid.hpp:1215
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition: CpGrid.hpp:1212
[ provides Dune::Grid ]
Definition: CpGrid.hpp:207
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.
Definition: CpGrid.hpp:801
std::string name() const
Get the grid name.
Definition: CpGrid.hpp:367
void switchToGlobalView()
Switch to the global view.
Definition: CpGrid.hpp:1478
const Vector & faceNormal(int face) const
Get the unit normal of a face.
Definition: CpGrid.hpp:1184
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:506
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:404
const InterfaceMap & pointScatterGatherInterface() const
Get an interface for gathering/scattering data attached to points with communication.
Definition: CpGrid.hpp:1472
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: CpGrid.hpp:603
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition: CpGrid.hpp:218
void gatherData(DataHandle &handle) const
Moves data from the distributed view to the global (all data on process) view.
Definition: CpGrid.hpp:1410
void globalRefine(int)
global refinement
Definition: CpGrid.hpp:517
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition: CpGrid.hpp:314
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: CpGrid.hpp:480
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:393
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:383
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition: CpGrid.hpp:1315
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGrid.hpp:354
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:424
void processEclipseFormat(const grdecl &input_data, bool remove_ij_boundary, bool turn_normals=false)
Read the Eclipse grid format ('grdecl').
Definition: CpGrid.cpp:551
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition: CpGrid.hpp:580
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition: CpGrid.hpp:1165
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: CpGrid.hpp:494
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: CpGrid.hpp:487
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition: CpGrid.hpp:1031
int maxLevel() const
Return maximum level defined in this grid.
Definition: CpGrid.hpp:375
const cpgrid::OrientedEntityTable< 0, 1 >::row_type cellFaceRow(int cell) const
Get a list of indices identifying all faces of a cell.
Definition: CpGrid.hpp:981
double faceArea(int face) const
Get the area of a face.
Definition: CpGrid.hpp:1171
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition: CpGrid.hpp:1177
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.
Definition: CpGrid.hpp:765
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.
Definition: CpGrid.hpp:659
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition: CpGrid.hpp:455
cpgrid::Entity< codim > entity(const cpgrid::EntityPointer< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition: CpGrid.hpp:529
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition: CpGrid.hpp:1043
bool loadBalance(int overlapLayers=1, bool useZoltan=true)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:631
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:448
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.
Definition: CpGrid.hpp:829
unsigned int ghostSize(int, int) const
Size of the ghost cell layer on a given level.
Definition: CpGrid.hpp:598
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.
Definition: CpGrid.hpp:692
int numVertices() const
Get The number of vertices.
Definition: CpGrid.hpp:956
void getIJK(const int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGrid.hpp:338
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition: CpGrid.hpp:974
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGrid.hpp:347
const CollectiveCommunication & comm() const
Get the collective communication object.
Definition: CpGrid.hpp:920
void createCartesian(const std::array< int, 3 > &dims, const std::array< double, 3 > &cellsize)
Create a cartesian grid.
Definition: CpGrid.cpp:448
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition: CpGrid.hpp:1196
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.hpp:995
int numCellFaces(int cell) const
Get the number of faces of a cell.
Definition: CpGrid.hpp:966
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:432
unsigned int overlapSize(int, int) const
Size of the overlap on a given level.
Definition: CpGrid.hpp:592
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:516
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition: CpGrid.hpp:586
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:414
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir) const
The new communication interface.
Definition: CpGrid.hpp:914
CpGrid()
Default constructor.
Definition: CpGrid.cpp:119
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition: CpGrid.hpp:1253
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data attached to cells with communication.
Definition: CpGrid.hpp:1465
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: CpGrid.hpp:501
void switchToDistributedView()
Switch to the distributed view.
Definition: CpGrid.hpp:1484
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition: CpGrid.hpp:1259
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.
Definition: CpGrid.hpp:882
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:440
int numCells() const
Get the number of cells.
Definition: CpGrid.hpp:946
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.
Definition: CpGrid.hpp:857
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: CpGrid.hpp:471
int numFaces() const
Get the number of faces.
Definition: CpGrid.hpp:951
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.
Definition: CpGrid.hpp:721
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition: CpGrid.hpp:1049
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition: CpGrid.hpp:901
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: CpGrid.hpp:510
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGrid.hpp:464
double cellVolume(int cell) const
Get the volume of the cell.
Definition: CpGrid.hpp:1190
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition: CpGrid.hpp:1390
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells.
Definition: CpGrid.hpp:326
std::map< int, std::list< int > > InterfaceMap
The type of the map describing communication interfaces.
Definition: CpGrid.hpp:1434
Struct that hods all the data needed to represent a Cpgrid.
Definition: CpGridData.hpp:123
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition: CpGridData.hpp:602
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGridData.hpp:252
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGridData.cpp:145
void getIJK(int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGridData.hpp:238
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGridData.hpp:259
Represents an entity of a given codim, with positive or negative orientation.
Definition: EntityRep.hpp:98
Definition: Entity.hpp:70
This class encapsulates geometry for both vertices, intersections and cells.
Definition: Geometry.hpp:68
The global id set for Dune.
Definition: Indexsets.hpp:325
Only needs to provide interface for doing nothing.
Definition: Iterators.hpp:89
Definition: Indexsets.hpp:54
Definition: Intersection.hpp:291
Definition: Intersection.hpp:68
Iterator intended to be used as LeafIterator and LevelIterator (no difference due to no adaptivity) f...
Definition: Iterators.hpp:56
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
int size() const
Returns the number of rows in the table.
Definition: SparseTable.hpp:121
int dataSize() const
Returns the number of data elements.
Definition: SparseTable.hpp:141
int size() const
Returns the number of rows in the table.
Definition: SparseTable.hpp:121
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
EdgeWeightMethod
enum for choosing Methods for weighting graph-edges correspoding to cell interfaces in Zoltan's graph...
Definition: GridEnums.hpp:34
@ defaultTransEdgeWgt
Use the transmissibilities as edge weights.
Definition: GridEnums.hpp:38
FieldVector< T, 3 > cross(const FieldVector< T, 3 > &a, const FieldVector< T, 3 > &b)
Definition: Volumes.hpp:58
Low-level corner-point processing routines and supporting data structures.
face_tag
Connection taxonomy.
Definition: preprocess.h:66
@ K_FACE
Connection topologically normal to I-J plane.
Definition: preprocess.h:69
@ J_FACE
Connection topologically normal to I-K plane.
Definition: preprocess.h:68
@ NNC_FACE
Arbitrary non-neighbouring connection.
Definition: preprocess.h:70
@ I_FACE
Connection topologically normal to J-K plane.
Definition: preprocess.h:67
Definition: CpGrid.hpp:194
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:147
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition.
Definition: CpGrid.hpp:149
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition: CpGrid.hpp:151
Traits associated with a specific codim.
Definition: CpGrid.hpp:120
cpgrid::Entity< cd > Entity
The type of the entity.
Definition: CpGrid.hpp:129
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition: CpGrid.hpp:123
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition: CpGrid.hpp:126
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition: CpGrid.hpp:135
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition: CpGrid.hpp:132
cpgrid::EntityPointer< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:141
cpgrid::EntityPointer< cd > EntityPointer
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:138
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:159
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:163
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:161
Definition: CpGrid.hpp:100
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition: CpGrid.hpp:173
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition: CpGrid.hpp:109
Dune::GridView< DefaultLeafGridViewTraits< CpGrid > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:170
Dune::GridView< DefaultLevelGridViewTraits< CpGrid > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:168
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition: CpGrid.hpp:177
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition: CpGrid.hpp:111
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition: CpGrid.hpp:175
GlobalIdSet LocalIdSet
The type of the local id set.
Definition: CpGrid.hpp:179
Dune::MPIHelper::MPICommunicator MPICommunicator
The type of the collective communication.
Definition: CpGrid.hpp:183
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition: CpGrid.hpp:114
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition: CpGrid.hpp:107
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition: CpGrid.hpp:105
CpGrid Grid
The type that implements the grid.
Definition: CpGrid.hpp:102
Raw corner-point specification of a particular geological model.
Definition: preprocess.h:56