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 #include <dune/common/parallel/variablesizecommunicator.hh>
56 #include <dune/grid/common/capabilities.hh>
57 #include <dune/grid/common/grid.hh>
58 #include <dune/grid/common/gridenums.hh>
60 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
62 #include "cpgrid/Intersection.hpp"
63 #include "cpgrid/Entity.hpp"
64 #include "cpgrid/Geometry.hpp"
65 #include "cpgrid/Iterators.hpp"
66 #include "cpgrid/Indexsets.hpp"
67 #include "cpgrid/DefaultGeometryPolicy.hpp"
68 #include "common/Volumes.hpp"
71 #include <opm/grid/utility/OpmParserIncludes.hpp>
137 template <PartitionIteratorType pitype>
149 template <PartitionIteratorType pitype>
152 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 5)
153 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> >
LevelGridView;
156 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> >
LeafGridView;
158 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid, pitype> >
LevelGridView;
161 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid, pitype> >
LeafGridView;
166 #if DUNE_VERSION_NEWER(DUNE_GRID, 2, 5)
167 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid> >
LevelGridView;
170 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid> >
LeafGridView;
172 typedef Dune::GridView<DefaultLevelGridViewTraits<CpGrid, Dune::All_Partition> >
LevelGridView;
175 typedef Dune::GridView<DefaultLeafGridViewTraits<CpGrid, Dune::All_Partition> >
LeafGridView;
190 typedef Dune::CollectiveCommunication<MPICommunicator> CollectiveCommunication;
214 :
public GridDefaultImplementation<3, 3, double, CpGridFamily >
218 typedef cpgrid::OpmEclipseStateType OpmEclipseStateType;
250 void processEclipseFormat(
const Opm::EclipseGrid& ecl_grid,
bool periodic_extension,
bool turn_normals =
false,
bool clip_z =
false,
259 const std::vector<double>& poreVolume = std::vector<double>());
281 const array<double, 3>& cellsize);
288 return current_view_data_->logical_cartesian_size_;
300 return current_view_data_->global_cell_;
310 void getIJK(
const int c, std::array<int,3>& ijk)
const
312 current_view_data_->
getIJK(c, ijk);
355 typename Traits::template Codim<codim>::LevelIterator
lbegin (
int level)
const
358 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
365 typename Traits::template Codim<codim>::LevelIterator
lend (
int level)
const
368 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
375 template<
int codim, PartitionIteratorType PiType>
376 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator
lbegin (
int level)
const
379 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
385 template<
int codim, PartitionIteratorType PiType>
386 typename Traits::template Codim<codim>::template Partition<PiType>::LevelIterator
lend (
int level)
const
389 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
396 typename Traits::template Codim<codim>::LeafIterator
leafbegin()
const
404 typename Traits::template Codim<codim>::LeafIterator
leafend()
const
411 template<
int codim, PartitionIteratorType PiType>
412 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator
leafbegin()
const
419 template<
int codim, PartitionIteratorType PiType>
420 typename Traits::template Codim<codim>::template Partition<PiType>::LeafIterator
leafend()
const
427 int size (
int level,
int codim)
const
430 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
438 return current_view_data_->
size(codim);
443 int size (
int level, GeometryType type)
const
446 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
452 int size (GeometryType type)
const
454 return current_view_data_->
size(type);
461 return *current_view_data_->global_id_set_;
468 return *current_view_data_->local_id_set_;
476 DUNE_THROW(GridError,
"levelIndexSet of nonexisting level " << level <<
" requested!");
477 return *current_view_data_->index_set_;
484 return *current_view_data_->index_set_;
491 std::cout <<
"Warning: Global refinement not implemented, yet." << std::endl;
494 const std::vector< Dune :: GeometryType >& geomTypes(
const int codim )
const
579 return current_view_data_->unique_boundary_ids_.size();
583 unsigned int numBndSegs = 0;
585 for (
int i = 0; i < num_faces; ++i) {
587 if (current_view_data_->face_to_cell_[face].
size() == 1) {
604 return get<0>(scatterGrid(
nullptr,
nullptr, overlapLayers ));
619 std::pair<bool, std::unordered_set<std::string> >
621 const double* transmissibilities =
nullptr,
624 return scatterGrid(ecl, transmissibilities, overlapLayers);
640 template<
class DataHandle>
641 std::pair<bool, std::unordered_set<std::string> >
643 const OpmEclipseStateType* ecl,
644 const double* transmissibilities =
nullptr,
647 auto ret =
loadBalance(ecl, transmissibilities, overlapLayers);
658 template<
class DataHandle>
674 template<
class DataHandle>
675 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir,
int )
const
687 template<
class DataHandle>
688 void communicate (DataHandle& data, InterfaceType iftype, CommunicationDirection dir)
const
690 current_view_data_->
communicate(data, iftype, dir);
694 const CollectiveCommunication&
comm ()
const
696 return current_view_data_->ccobj_;
710 typedef Dune::FieldVector<double, 3> Vector;
713 const std::vector<double>& zcornData()
const {
714 return data_->zcornData();
722 return current_view_data_->cell_to_face_.
size();
727 return current_view_data_->face_to_cell_.
size();
732 return current_view_data_->geomVector<3>().
size();
750 return current_view_data_->cell_to_face_[
cpgrid::EntityRep<0>(cell,
true)][local_index].index();
776 bool a = (local_index == 0);
777 bool b = r[0].orientation();
778 bool use_first = a ? b : !b;
780 int r_size = r.size();
783 if(r[0].index()==std::numeric_limits<int>::max()){
788 if(r.size()>1 && r[1].index()==std::numeric_limits<int>::max())
794 return use_first ? r[0].index() : r[1].index();
796 return use_first ? r[index].index() : -1;
807 return current_view_data_->cell_to_face_.
dataSize();
809 int numFaceVertices(
int face)
const
811 return current_view_data_->face_to_point_[face].
size();
819 return current_view_data_->face_to_point_[face][local_index];
828 const int nv = current_view_data_->cell_to_point_[cell_index].size();
830 for (
int i=0; i<nv; ++i) {
831 zz +=
vertexPosition(current_view_data_->cell_to_point_[cell_index][i])[nd-1];
836 const Vector faceCenterEcl(
int cell_index,
int faceTag)
const
851 static const int faceVxMap[ 6 ][ 4 ] = { {0, 2, 4, 6},
860 assert (current_view_data_->cell_to_point_[cell_index].size() == 8);
862 for(
int i=0; i<4; ++i )
864 center +=
vertexPosition(current_view_data_->cell_to_point_[cell_index][ faceVxMap[ faceTag ][ i ] ]);
867 for (
int i=0; i<3; ++i) {
874 const Vector faceAreaNormalEcl(
int face)
const
877 const int nd = Vector::dimension;
878 const int nv = numFaceVertices(face);
892 Vector areaNormal =
cross(a,b);
893 for (
int i=0; i<nd; ++i) {
903 Vector areaNormal =
cross(a,b);
911 int k = (nv % 2) ? 0 : nv - 1;
913 Vector areaNormal(0.0);
915 for (
int i = 1; i < h; ++i)
918 Vector b =
vertexPosition(current_view_data_->face_to_point_[face][2*i+1]) -
vertexPosition(current_view_data_->face_to_point_[face][2*i-1]);
919 areaNormal +=
cross(a,b);
924 Vector b =
vertexPosition(current_view_data_->face_to_point_[face][k]) -
vertexPosition(current_view_data_->face_to_point_[face][2*h-1]);
925 areaNormal +=
cross(a,b);
960 return current_view_data_->face_normals_.get(face);
979 :
public RandomAccessIteratorFacade<CentroidIterator<codim>,
980 FieldVector<double, 3>,
981 const FieldVector<double, 3>&, int>
993 const FieldVector<double, 3>& dereference()
const
995 return iter_->center();
1001 const FieldVector<double, 3>& elementAt(
int n)
1003 return iter_[n]->center();
1039 int boundaryId(
int face)
const
1048 if (current_view_data_->face_to_cell_[f].
size() == 1) {
1051 ret = current_view_data_->unique_boundary_ids_[f];
1054 const bool normal_is_in =
1055 !(current_view_data_->face_to_cell_[f][0].orientation());
1056 enum face_tag tag = current_view_data_->face_tag_[f];
1060 ret = normal_is_in ? 1 : 2;
1064 ret = normal_is_in ? 3 : 4;
1069 ret = normal_is_in ? 5 : 6;
1083 template<
class Cell2FacesRowIterator>
1085 faceTag(
const Cell2FacesRowIterator& cell_face)
const
1096 const int cell = cell_face.getCellIndex();
1097 const int face = *cell_face;
1098 assert (0 <= cell); assert (cell <
numCells());
1099 assert (0 <= face); assert (face <
numFaces());
1104 const F2C& f2c = current_view_data_->face_to_cell_[f];
1105 const face_tag tag = current_view_data_->face_tag_[f];
1107 assert ((f2c.size() == 1) || (f2c.size() == 2));
1109 int inside_cell = 0;
1111 if ( f2c.size() == 2 )
1113 if ( f2c[1].index() == cell )
1118 const bool normal_is_in = ! f2c[inside_cell].orientation();
1123 return normal_is_in ? 0 : 1;
1127 return normal_is_in ? 2 : 3;
1132 return normal_is_in ? 4 : 5;
1134 OPM_THROW(std::logic_error,
"Unhandled face tag. This should never happen!");
1149 template<
class DataHandle>
1162 if(!distributed_data_)
1163 OPM_THROW(std::runtime_error,
"Moving Data only allowed with a load balanced grid!");
1164 distributed_data_->scatterData(handle, data_.get(), distributed_data_.get());
1177 template<
class DataHandle>
1181 if(!distributed_data_)
1182 OPM_THROW(std::runtime_error,
"Moving Data only allowed with a load balance grid!");
1183 distributed_data_->gatherData(handle, data_.get(), distributed_data_.get());
1190 typedef VariableSizeCommunicator<>::InterfaceMap
InterfaceMap;
1230 return *cell_scatter_gather_interfaces_;
1236 current_view_data_=data_.get();
1242 current_view_data_=distributed_data_.get();
1247 typedef cpgrid::CpGridData::ParallelIndexSet ParallelIndexSet;
1250 typedef cpgrid::CpGridData::RemoteIndices RemoteIndices;
1252 ParallelIndexSet& getCellIndexSet()
1254 return current_view_data_->cell_indexset_;
1257 RemoteIndices& getCellRemoteIndices()
1259 return current_view_data_->cell_remote_indices_;
1262 const ParallelIndexSet& getCellIndexSet()
const
1264 return current_view_data_->cell_indexset_;
1267 const RemoteIndices& getCellRemoteIndices()
const
1269 return current_view_data_->cell_remote_indices_;
1283 std::pair<bool, std::unordered_set<std::string> >
1284 scatterGrid(
const OpmEclipseStateType* ecl,
const double* transmissibilities,
1291 std::shared_ptr<cpgrid::CpGridData> data_;
1293 cpgrid::CpGridData* current_view_data_;
1295 std::shared_ptr<cpgrid::CpGridData> distributed_data_;
1301 std::shared_ptr<InterfaceMap> cell_scatter_gather_interfaces_;
1306 namespace Capabilities
1312 static const bool v =
true;
1319 static const bool v =
true;
1322 #if ! DUNE_VERSION_NEWER(DUNE_GRID, 2, 5)
1327 static const bool v =
true;
1334 static const bool v =
true;
1340 static const bool v =
true;
1347 static const bool v =
false;
1354 #include <dune/grid/cpgrid/PersistentContainer.hpp>
1355 #include <dune/grid/cpgrid/CartesianIndexMapper.hpp>
1356 #endif // OPM_CPGRID_HEADER
cpgrid::IdSet LocalIdSet
The type of the local id set.
Definition: CpGrid.hpp:185
Definition: CpGrid.hpp:201
int faceTag(const Cell2FacesRowIterator &cell_face) const
Get the cartesian tag associated with a face tag.
Definition: CpGrid.hpp:1085
Definition: CpGridData.hpp:93
Dune::MPIHelper::MPICommunicator MPICommunicator
The type of the collective communication.
Definition: CpGrid.hpp:189
void processEclipseFormat(const grdecl &input_data, double z_tolerance, bool remove_ij_boundary, bool turn_normals=false)
Read the Eclipse grid format ('grdecl').
Definition: CpGrid.cpp:276
cpgrid::Entity< codim > entity(const cpgrid::EntityPointer< codim > &seed) const
given an EntitySeed (or EntityPointer) return an entity object
Definition: CpGrid.hpp:501
std::vector< cpgrid::Geometry< 3-codim, 3 > >::const_iterator GeometryIterator
The type of the iterator over the codim geometries.
Definition: CpGrid.hpp:986
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:138
void switchToGlobalView()
Switch to the global view.
Definition: CpGrid.hpp:1234
Dune::GridView< DefaultLevelGridViewTraits< CpGrid, Dune::All_Partition > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:173
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.hpp:769
std::string name() const
Get the grid name.
Definition: CpGrid.hpp:339
cpgrid::HierarchicIterator HierarchicIterator
The type of the hierarchic iterator.
Definition: CpGrid.hpp:106
Raw corner-point specification of a particular geological model.
Definition: preprocess.h:56
cpgrid::IndexSet LeafIndexSet
The type of the leaf index set.
Definition: CpGrid.hpp:181
int numCellFaces(int cell) const
Get the number of faces of a cell.
Definition: CpGrid.hpp:740
cpgrid::GlobalIdSet GlobalIdSet
The type of the global id set.
Definition: CpGrid.hpp:183
Definition: CpGrid.hpp:91
Dune::GridView< DefaultLeafGridViewTraits< CpGrid, Dune::All_Partition > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:175
int cellFace(int cell, int local_index) const
Get a specific face of a cell.
Definition: CpGrid.hpp:748
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:755
std::pair< bool, std::unordered_set< std::string > > loadBalance(const OpmEclipseStateType *ecl, const double *transmissibilities=nullptr, int overlapLayers=1)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:620
void createCartesian(const array< int, 3 > &dims, const array< double, 3 > &cellsize)
Create a cartesian grid.
Definition: CpGrid.cpp:212
int numCellFaces() const
Get the sum of all faces attached to all cells.
Definition: CpGrid.hpp:805
void readSintefLegacyFormat(const std::string &grid_prefix)
Read the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:254
int numCells() const
Get the number of cells.
Definition: CpGrid.hpp:720
Definition: Indexsets.hpp:192
Definition: Intersection.hpp:67
cpgrid::Intersection LeafIntersection
The type of the intersection at the leafs of the grid.
Definition: CpGrid.hpp:97
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: CpGrid.hpp:443
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGrid.hpp:326
const std::vector< int > & globalCell() const
Retrieve mapping from internal ("compressed") active grid cells to external ("uncompressed") cells...
Definition: CpGrid.hpp:298
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:420
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGridData.hpp:208
double cellVolume(int cell) const
Get the volume of the cell.
Definition: CpGrid.hpp:964
[ provides Dune::Grid ]
Definition: CpGrid.hpp:213
unsigned int ghostSize(int, int) const
Size of the ghost cell layer on a given level.
Definition: CpGrid.hpp:570
Traits associated with a specific grid partition type.
Definition: CpGrid.hpp:150
cpgrid::IntersectionIterator LevelIntersectionIterator
The type of the intersection iterator at the levels of the grid.
Definition: CpGrid.hpp:103
This class encapsulates geometry for both vertices, intersections and cells.
Definition: CpGridData.hpp:92
Traits::template Codim< codim >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:355
Iterator intended to be used as LeafIterator and LevelIterator (no difference due to no adaptivity) f...
Definition: Entity.hpp:52
int numFaces() const
Get the number of faces.
Definition: CpGrid.hpp:725
unsigned int ghostSize(int) const
Size of the ghost cell layer on the leaf level.
Definition: CpGrid.hpp:558
std::map< int, std::list< int > > InterfaceMap
The type of the map describing communication interfaces.
Definition: CpGrid.hpp:1197
const InterfaceMap & cellScatterGatherInterface() const
Get an interface for gathering/scattering data with communication.
Definition: CpGrid.hpp:1228
unsigned int numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: CpGrid.hpp:575
void setUniqueBoundaryIds(bool uids)
Set whether we want to have unique boundary ids.
Definition: CpGridData.hpp:215
cpgrid::Entity< cd > Entity
The type of the entity.
Definition: CpGrid.hpp:121
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir)
communicate objects for all codims on a given level
Definition: CpGridData.hpp:512
bool loadBalance(DataHandle &data, int overlapLayers=1)
Distributes this grid and data over the available nodes in a distributed machine. ...
Definition: CpGrid.hpp:659
Traits::template Codim< codim >::LeafIterator leafend() const
one past the end of the sequence of leaf entities
Definition: CpGrid.hpp:404
CpGrid Grid
The type that implements the grid.
Definition: CpGrid.hpp:94
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGrid.hpp:436
int maxLevel() const
Return maximum level defined in this grid.
Definition: CpGrid.hpp:347
cpgrid::Iterator< cd, pitype > LevelIterator
The type of the iterator over the level entities of this codim on this partition. ...
Definition: CpGrid.hpp:141
const CollectiveCommunication & comm() const
Get the collective communication object.
Definition: CpGrid.hpp:694
Traits::template Codim< codim >::template Partition< PiType >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:412
void switchToDistributedView()
Switch to the distributed view.
Definition: CpGrid.hpp:1240
cpgrid::IntersectionIterator LeafIntersectionIterator
The type of the intersection iterator at the leafs of the grid.
Definition: CpGrid.hpp:101
CentroidIterator< 1 > beginFaceCentroids() const
Get an iterator over the face centroids positioned at the first one.
Definition: CpGrid.hpp:1033
CpGrid()
Default constructor.
Definition: CpGrid.cpp:58
Traits associated with a specific codim.
Definition: CpGrid.hpp:111
const Traits::GlobalIdSet & globalIdSet() const
Access to the GlobalIdSet.
Definition: CpGrid.hpp:459
CentroidIterator< 0 > beginCellCentroids() const
Get an iterator over the cell centroids positioned at the first one.
Definition: CpGrid.hpp:1027
int numVertices() const
Get The number of vertices.
Definition: CpGrid.hpp:730
unsigned int overlapSize(int) const
Size of the overlap on the leaf level.
Definition: CpGrid.hpp:552
cpgrid::Iterator< cd, All_Partition > LeafIterator
The type of the iterator over all leaf entities of this codim.
Definition: CpGrid.hpp:127
const Vector & faceCentroid(int face) const
Get the coordinates of the center of a face.
Definition: CpGrid.hpp:951
double faceArea(int face) const
Get the area of a face.
Definition: CpGrid.hpp:945
Traits::template Codim< codim >::LeafIterator leafbegin() const
Iterator to first leaf entity of given codim.
Definition: CpGrid.hpp:396
int dataSize() const
Returns the number of data elements.
Definition: SparseTable.hpp:144
Connection topologically parallel to I-K plane.
Definition: preprocess.h:69
int size(int codim) const
number of leaf entities per codim in this process
Definition: CpGridData.cpp:132
const Traits::LocalIdSet & localIdSet() const
Access to the LocalIdSet.
Definition: CpGrid.hpp:466
An iterator over the centroids of the geometry of the entities.
Definition: CpGrid.hpp:978
std::pair< bool, std::unordered_set< std::string > > loadBalance(DataHandle &data, const OpmEclipseStateType *ecl, const double *transmissibilities=nullptr, int overlapLayers=1)
Distributes this grid and data over the available nodes in a distributed machine. ...
Definition: CpGrid.hpp:642
const Traits::LevelIndexSet & levelIndexSet(int level) const
Access to the LevelIndexSets.
Definition: CpGrid.hpp:473
Only needs to provide interface for doing nothing.
Definition: Iterators.hpp:88
cpgrid::EntityPointer< cd > EntitySeed
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:133
Definition: Indexsets.hpp:250
int size() const
Returns the number of rows in the table.
Definition: SparseTable.hpp:124
unsigned int overlapSize(int, int) const
Size of the overlap on a given level.
Definition: CpGrid.hpp:564
int size(int level, int codim) const
Number of grid entities per level and codim.
Definition: CpGrid.hpp:427
Connection topologically parallel to I-J plane.
Definition: preprocess.h:70
cpgrid::Iterator< cd, pitype > LeafIterator
The type of the iterator over the leaf entities of this codim on this partition.
Definition: CpGrid.hpp:143
const std::array< int, 3 > & logicalCartesianSize() const
The logical cartesian size of the global grid.
Definition: CpGrid.hpp:286
cpgrid::Geometry< 3-cd, 3 > LocalGeometry
The type of the local geometry associated with the entity.
Definition: CpGrid.hpp:118
void getIJK(int c, std::array< int, 3 > &ijk) const
Extract Cartesian index triplet (i,j,k) of an active cell.
Definition: CpGridData.hpp:194
face_tag
Connection taxonomy.
Definition: preprocess.h:67
cpgrid::EntityPointer< cd > EntityPointer
The type of the entity pointer for entities of this codim.
Definition: CpGrid.hpp:130
CpGridFamily GridFamily
Family typedef, why is this not defined by Grid<>?
Definition: CpGrid.hpp:226
Dune::GridView< DefaultLeafGridViewTraits< CpGrid, pitype > > LeafGridView
The type of the leaf grid view associated with this partition type.
Definition: CpGrid.hpp:161
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write the Sintef legacy grid format ('topogeom').
Definition: CpGrid.cpp:258
cpgrid::IndexSet LevelIndexSet
The type of the level index set.
Definition: CpGrid.hpp:179
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: CpGrid.hpp:452
CentroidIterator(GeometryIterator iter)
Constructs a new iterator from an iterator over the geometries.
Definition: CpGrid.hpp:989
bool uniqueBoundaryIds() const
Is the grid currently using unique boundary ids?
Definition: CpGrid.hpp:319
void globalRefine(int)
global refinement
Definition: CpGrid.hpp:489
cpgrid::Geometry< 3-cd, 3 > Geometry
The type of the geometry associated with the entity.
Definition: CpGrid.hpp:115
const Traits::LeafIndexSet & leafIndexSet() const
Access to the LeafIndexSet.
Definition: CpGrid.hpp:482
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int) const
The new communication interface.
Definition: CpGrid.hpp:675
Low-level corner-point processing routines and supporting data structures.
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:310
Represents an entity of a given codim, with positive or negative orientation.
Definition: CpGridData.hpp:94
const Vector & vertexPosition(int vertex) const
Get the Position of a vertex.
Definition: CpGrid.hpp:939
Definition: Indexsets.hpp:52
cpgrid::Intersection LevelIntersection
The type of the intersection at the levels of the grid.
Definition: CpGrid.hpp:99
Definition: Intersection.hpp:286
void scatterData(DataHandle &handle)
Moves data from the global (all data on process) view to the distributed view.
Definition: CpGrid.hpp:1159
void communicate(DataHandle &data, InterfaceType iftype, CommunicationDirection dir) const
The new communication interface.
Definition: CpGrid.hpp:688
int faceVertex(int face, int local_index) const
Get the index identifying a vertex of a face.
Definition: CpGrid.hpp:817
Traits::template Codim< codim >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:365
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
void gatherData(DataHandle &handle)
Moves data from the distributed view to the global (all data on process) view.
Definition: CpGrid.hpp:1178
const Vector & faceNormal(int face) const
Get the unit normal of a face.
Definition: CpGrid.hpp:958
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lend(int level) const
one past the end on this level
Definition: CpGrid.hpp:386
Connection topologically parallel to J-K plane.
Definition: preprocess.h:68
FieldVector< T, 3 > cross(const FieldVector< T, 3 > &a, const FieldVector< T, 3 > &b)
Definition: Volumes.hpp:58
Traits::template Codim< codim >::template Partition< PiType >::LevelIterator lbegin(int level) const
Iterator to first entity of given codim on level.
Definition: CpGrid.hpp:376
Dune::GridView< DefaultLevelGridViewTraits< CpGrid, pitype > > LevelGridView
The type of the level grid view associated with this partition type.
Definition: CpGrid.hpp:159
cpgrid::Iterator< cd, All_Partition > LevelIterator
The type of the iterator over all level entities of this codim.
Definition: CpGrid.hpp:124
bool loadBalance(int overlapLayers=1)
Distributes this grid over the available nodes in a distributed machine.
Definition: CpGrid.hpp:601
double cellCenterDepth(int cell_index) const
Get vertical position of cell center ("zcorn" average).
Definition: CpGrid.hpp:823
const Vector & cellCentroid(int cell) const
Get the coordinates of the center of a cell.
Definition: CpGrid.hpp:970