3#ifndef DUNE_POLYHEDRALGRID_GRID_HH
4#define DUNE_POLYHEDRALGRID_GRID_HH
10#include <opm/grid/utility/platform_dependent/disable_warnings.h>
13#include <dune/common/version.hh>
14#include <dune/common/parallel/mpihelper.hh>
17#include <dune/grid/common/grid.hh>
19#include <dune/common/parallel/communication.hh>
22#include <opm/grid/polyhedralgrid/capabilities.hh>
23#include <opm/grid/polyhedralgrid/declaration.hh>
24#include <opm/grid/polyhedralgrid/entity.hh>
25#include <opm/grid/polyhedralgrid/entityseed.hh>
26#include <opm/grid/polyhedralgrid/geometry.hh>
27#include <opm/grid/polyhedralgrid/gridview.hh>
28#include <opm/grid/polyhedralgrid/idset.hh>
31#include <opm/grid/utility/platform_dependent/reenable_warnings.h>
33#include <opm/common/ErrorMacros.hpp>
38#include <opm/grid/GridManager.hpp>
40#include <opm/grid/MinpvProcessor.hpp>
43#include <opm/input/eclipse/EclipseState/Grid/EclipseGrid.hpp>
53 template<
int dim,
int dimworld,
typename coord_t >
60 typedef coord_t ctype;
68 static const int dimension = dim;
69 static const int dimensionworld = dimworld;
71 typedef Dune::FieldVector< ctype, dimensionworld > GlobalCoordinate ;
78 typedef Dune::Intersection< const Grid, LeafIntersectionImpl > LeafIntersection;
79 typedef Dune::Intersection< const Grid, LevelIntersectionImpl > LevelIntersection;
81 typedef Dune::IntersectionIterator< const Grid, LeafIntersectionIteratorImpl, LeafIntersectionImpl > LeafIntersectionIterator;
82 typedef Dune::IntersectionIterator< const Grid, LevelIntersectionIteratorImpl, LevelIntersectionImpl > LevelIntersectionIterator;
85 typedef Dune::EntityIterator< 0, const Grid, HierarchicIteratorImpl > HierarchicIterator;
97 typedef Dune::Entity< codim, dimension, const Grid, PolyhedralGridEntity > Entity;
100 typedef Entity EntityPointer;
105 template< PartitionIteratorType pitype >
109 typedef Dune::EntityIterator< codim, const Grid, LeafIteratorImpl > LeafIterator;
111 typedef LeafIterator LevelIterator;
114 typedef typename Partition< All_Partition >::LeafIterator LeafIterator;
115 typedef typename Partition< All_Partition >::LevelIterator LevelIterator;
124 typedef Dune::MPIHelper::MPICommunicator MPICommunicator;
125 using Communication = Dune::Communication<MPICommunicator>;
126 using CollectiveCommunication = Dune::Communication<MPICommunicator>;
128 template< PartitionIteratorType pitype >
131 typedef Dune::GridView< PolyhedralGridViewTraits< dim, dimworld, ctype, pitype > > LeafGridView;
132 typedef Dune::GridView< PolyhedralGridViewTraits< dim, dimworld, ctype, pitype > > LevelGridView;
135 typedef typename Partition< All_Partition >::LeafGridView LeafGridView;
136 typedef typename Partition< All_Partition >::LevelGridView LevelGridView;
153 template <
int dim,
int dimworld,
typename coord_t >
156 :
public GridDefaultImplementation
157 < dim, dimworld, coord_t, PolyhedralGridFamily< dim, dimworld, coord_t > >
162 typedef GridDefaultImplementation
178 typedef std::unique_ptr< UnstructuredGridType, UnstructuredGridDeleter > UnstructuredGridPtr;
180 static UnstructuredGridPtr
181 allocateGrid ( std::size_t nCells, std::size_t nFaces, std::size_t nFaceNodes, std::size_t nCellFaces, std::size_t nNodes )
186 DUNE_THROW( GridError,
"Unable to allocate grid" );
187 return UnstructuredGridPtr( grid );
191 computeGeometry ( UnstructuredGridPtr& ug )
201 typedef PolyhedralGridFamily< dim, dimworld, coord_t > GridFamily;
208 typedef typename GridFamily::Traits
Traits;
216 template<
int codim >
237 template< PartitionIteratorType pitype >
240 typedef typename GridFamily::Traits::template Partition< pitype >::LevelGridView LevelGridView;
241 typedef typename GridFamily::Traits::template Partition< pitype >::LeafGridView LeafGridView;
246 typedef typename Partition< All_Partition >::LeafGridView LeafGridView;
308 typedef typename Traits::ctype
ctype;
314 typedef typename Traits :: GlobalCoordinate GlobalCoordinate;
328 const std::vector<double>& poreVolumes = std::vector<double> ())
329 : gridPtr_( createGrid( inputGrid, poreVolumes ) ),
331 comm_( MPIHelper::getCommunicator() ),
332 leafIndexSet_( *this ),
333 globalIdSet_( *this ),
334 localIdSet_( *this ),
347 const std::vector< double >& dx )
348 : gridPtr_( createGrid( n, dx ) ),
350 comm_( MPIHelper::getCommunicator()),
351 leafIndexSet_( *this ),
352 globalIdSet_( *this ),
353 localIdSet_( *this ),
366 : gridPtr_( std::move( gridPtr ) ),
368 comm_( MPIHelper::getCommunicator() ),
369 leafIndexSet_( *this ),
370 globalIdSet_( *this ),
371 localIdSet_( *this ),
387 comm_( MPIHelper::getCommunicator() ),
388 leafIndexSet_( *this ),
389 globalIdSet_( *this ),
390 localIdSet_( *this ),
400 operator const UnstructuredGridType& ()
const {
return grid_; }
427 int size (
int ,
int codim )
const
429 return size( codim );
444 else if ( codim == 1 )
448 else if ( codim == dim )
454 std::cerr <<
"Warning: codimension " << codim <<
" not available in PolyhedralGrid" << std::endl;
467 int size (
int , GeometryType type )
const
469 return size( dim - type.dim() );
476 int size ( GeometryType type )
const
478 return size( dim - type.dim() );
489 return nBndSegments_;
493 template<
int codim >
496 return leafbegin< codim, All_Partition >();
499 template<
int codim >
502 return leafend< codim, All_Partition >();
505 template<
int codim, PartitionIteratorType pitype >
506 typename Codim< codim >::template Partition< pitype >::LeafIterator
509 typedef typename Traits::template Codim< codim >::template Partition< pitype >::LeafIteratorImpl Impl;
510 return Impl( extraData(),
true );
513 template<
int codim, PartitionIteratorType pitype >
514 typename Codim< codim >::template Partition< pitype >::LeafIterator
517 typedef typename Traits::template Codim< codim >::template Partition< pitype >::LeafIteratorImpl Impl;
518 return Impl( extraData(),
false );
521 template<
int codim >
524 return leafbegin< codim, All_Partition >();
527 template<
int codim >
530 return leafend< codim, All_Partition >();
533 template<
int codim, PartitionIteratorType pitype >
534 typename Codim< codim >::template Partition< pitype >::LevelIterator
535 lbegin (
const int )
const
537 return leafbegin< codim, pitype > ();
540 template<
int codim, PartitionIteratorType pitype >
541 typename Codim< codim >::template Partition< pitype >::LevelIterator
542 lend (
const int )
const
544 return leafend< codim, pitype > ();
559 return leafIndexSet();
564 return leafIndexSet_;
567 void globalRefine (
int )
597 template<
class DataHandle >
626 return (codim == 0 ) ? 1 : 0;
662 template<
class DataHandle>
665 CommunicationDirection ,
668 OPM_THROW(std::runtime_error,
"communicate not implemented for polyhedreal grid!");
683 template<
class DataHandle>
686 CommunicationDirection )
const
688 OPM_THROW(std::runtime_error,
"communicate not implemented for polyhedreal grid!");
694 OPM_THROW(std::runtime_error,
"switch to global view not implemented for polyhedreal grid!");
700 OPM_THROW(std::runtime_error,
"switch to distributed view not implemented for polyhedreal grid!");
711 const CommunicationType &
comm ()
const
747 template<
class DataHandle,
class Data >
767 template<
class DofManager >
774 template< PartitionIteratorType pitype >
777 typedef typename Partition< pitype >::LevelGridView View;
778 typedef typename View::GridViewImp ViewImp;
779 return View( ViewImp( *
this ) );
783 template< PartitionIteratorType pitype >
786 typedef typename Traits::template Partition< pitype >::LeafGridView View;
787 typedef typename View::GridViewImp ViewImp;
788 return View( ViewImp( *
this ) );
794 typedef typename LevelGridView::GridViewImp ViewImp;
801 typedef typename LeafGridView::GridViewImp ViewImp;
802 return LeafGridView( ViewImp( *
this ) );
806 template<
class EntitySeed >
813 return EntityPointer( EntityPointerImpl( EntityImpl( extraData(), seed ) ) );
817 template<
class EntitySeed >
822 return EntityImpl( extraData(), seed );
844 const std::array<int, 3>& logicalCartesianSize()
const
849 const int* globalCell()
const
855 const int* globalCellPtr()
const
860 void getIJK(
const int c, std::array<int,3>& ijk)
const
862 int gc = globalCell()[c];
863 ijk[0] = gc % logicalCartesianSize()[0]; gc /= logicalCartesianSize()[0];
864 ijk[1] = gc % logicalCartesianSize()[1];
865 ijk[2] = gc / logicalCartesianSize()[1];
882 template<
class DataHandle>
885 OPM_THROW(std::runtime_error,
"ScatterData not implemented for polyhedral grid!");
890 UnstructuredGridType* createGrid(
const Opm::EclipseGrid& inputGrid,
const std::vector< double >& poreVolumes )
const
894 g.
dims[0] = inputGrid.getNX();
895 g.dims[1] = inputGrid.getNY();
896 g.dims[2] = inputGrid.getNZ();
898 std::vector<double>
coord = inputGrid.getCOORD( );
899 std::vector<double>
zcorn = inputGrid.getZCORN( );
900 std::vector<int>
actnum = inputGrid.getACTNUM( );
902 g.coord =
coord.data();
903 g.zcorn =
zcorn.data();
906 if (!poreVolumes.empty() && (inputGrid.getMinpvMode() != Opm::MinpvMode::Inactive))
909 const std::vector<double>& minpvv = inputGrid.getMinpvVector();
914 const size_t cartGridSize = g.dims[0] * g.dims[1] * g.dims[2];
915 std::vector<double> thickness(cartGridSize);
916 for (
size_t i = 0; i < cartGridSize; ++i) {
917 thickness[i] = inputGrid.getCellThickness(i);
919 const double z_tolerance = inputGrid.isPinchActive() ? inputGrid.getPinchThresholdThickness() : 0.0;
920 mp.process(thickness, z_tolerance, poreVolumes, minpvv,
actnum, opmfil,
zcorn.data());
934 const double z_tolerance = inputGrid.isPinchActive() ?
935 inputGrid.getPinchThresholdThickness() : 0.0;
938 OPM_THROW(std::runtime_error,
"Failed to construct grid.");
944 UnstructuredGridType* createGrid(
const std::vector< int >& n,
const std::vector< double >& dx )
const
946 UnstructuredGridType* cgrid = nullptr ;
947 assert(
int(n.size()) == dim );
959 OPM_THROW(std::runtime_error,
"Failed to construct grid.");
965#if DUNE_VERSION_LT_REV(DUNE_GRID, 2, 7, 1)
966 using Base::getRealImplementation;
968 typedef typename Traits :: ExtraData ExtraData;
969 ExtraData extraData ()
const {
return this; }
971 template <
class EntitySeed>
972 int corners(
const EntitySeed& seed )
const
974 const int codim = EntitySeed :: codimension;
975 const int index = seed.index();
977 return cellVertices_[ index ].size();
985 template <
class EntitySeed>
987 corner (
const EntitySeed& seed,
const int i )
const
989 const int codim = EntitySeed :: codimension;
992 const int coordIndex = GlobalCoordinate :: dimension * cellVertices_[ seed.index() ][ i ];
1000 const int crners = corners( seed );
1001 const int crner = (crners == 4 && EntitySeed :: dimension == 3 && i > 1 ) ? 5 - i : i;
1003 return copyToGlobalCoordinate( grid_.
node_coordinates + GlobalCoordinate :: dimension * faceVertex );
1007 const int coordIndex = GlobalCoordinate :: dimension * seed.index();
1010 return GlobalCoordinate( 0 );
1013 template <
class EntitySeed>
1014 int subEntities(
const EntitySeed& seed,
const int codim )
const
1016 const int index = seed.index();
1017 if( seed.codimension == 0 )
1024 return cellVertices_[ index ].size();
1026 else if( seed.codimension == 1 )
1033 else if ( seed.codimension == dim )
1041 template <
int codim,
class EntitySeedArg >
1042 typename Codim<codim>::EntitySeed
1043 subEntitySeed(
const EntitySeedArg& baseSeed,
const int i )
const
1045 assert( codim >= EntitySeedArg::codimension );
1046 assert( i>= 0 && i<subEntities( baseSeed, codim ) );
1047 typedef typename Codim<codim>::EntitySeed EntitySeed;
1050 if( codim == EntitySeedArg::codimension )
1052 return EntitySeed( baseSeed.index() );
1055 if( EntitySeedArg::codimension == 0 )
1061 else if ( codim == dim )
1063 return EntitySeed( cellVertices_[ baseSeed.index() ][ i ] );
1066 else if ( EntitySeedArg::codimension == 1 && codim == dim )
1071 DUNE_THROW(NotImplemented,
"codimension not available");
1072 return EntitySeed();
1075 template <
int codim>
1076 typename Codim<codim>::EntitySeed
1077 subEntitySeed(
const typename Codim<1>::EntitySeed& faceSeed,
const int i )
const
1079 assert( i>= 0 && i<subEntities( faceSeed, codim ) );
1080 typedef typename Codim<codim>::EntitySeed EntitySeed;
1083 return EntitySeed( faceSeed.index() );
1085 else if ( codim == dim )
1091 DUNE_THROW(NotImplemented,
"codimension not available");
1095 bool hasBoundaryIntersections(
const typename Codim<0>::EntitySeed& seed )
const
1097 const int faces = subEntities( seed, 1 );
1098 for(
int f=0; f<faces; ++f )
1100 const auto faceSeed = this->
template subEntitySeed<1>( seed, f );
1101 if( isBoundaryFace( faceSeed ) )
1107 bool isBoundaryFace(
const int face )
const
1110 const int facePos = 2 * face;
1114 bool isBoundaryFace(
const typename Codim<1>::EntitySeed& faceSeed )
const
1116 assert( faceSeed.isValid() );
1117 return isBoundaryFace( faceSeed.index() );
1120 int boundarySegmentIndex(
const typename Codim<0>::EntitySeed& seed,
const int face )
const
1122 const auto faceSeed = this->
template subEntitySeed<1>( seed, face );
1123 assert( faceSeed.isValid() );
1124 const int facePos = 2 * faceSeed.index();
1131 const std::vector< GeometryType > &geomTypes (
const unsigned int codim )
const
1133 static std::vector< GeometryType > emptyDummy;
1134 if (codim < geomTypes_.size())
1136 return geomTypes_[codim];
1142 template <
class Seed >
1143 GeometryType geometryType(
const Seed& seed )
const
1145 if( Seed::codimension == 0 )
1147 assert(!geomTypes( Seed::codimension ).empty());
1148 return geomTypes( Seed::codimension )[ 0 ];
1153 if( dim == 3 && Seed::codimension == 1 )
1156 const int nVx = corners( seed );
1158 face = Dune::GeometryTypes::cube(2);
1160 face = Dune::GeometryTypes::simplex(2);
1162 face = Dune::GeometryTypes::none(2);
1167 assert(!geomTypes( Seed::codimension ).empty());
1168 return geomTypes( Seed::codimension )[ 0 ];
1172 int indexInInside(
const typename Codim<0>::EntitySeed& seed,
const int i )
const
1174 return ( grid_.
cell_facetag ) ? cartesianIndexInInside( seed, i ) : i;
1177 int cartesianIndexInInside(
const typename Codim<0>::EntitySeed& seed,
const int i )
const
1179 assert( i>= 0 && i<subEntities( seed, 1 ) );
1183 typename Codim<0>::EntitySeed
1184 neighbor(
const typename Codim<0>::EntitySeed& seed,
const int i )
const
1186 const int face = this->
template subEntitySeed<1>( seed, i ).index();
1188 if( nb == seed.index() )
1193 typedef typename Codim<0>::EntitySeed EntitySeed;
1194 return EntitySeed( nb );
1198 indexInOutside(
const typename Codim<0>::EntitySeed& seed,
const int i )
const
1203 const int in_inside = cartesianIndexInInside( seed, i );
1204 return in_inside + ((in_inside % 2) ? -1 : 1);
1208 typedef typename Codim<0>::EntitySeed EntitySeed;
1209 EntitySeed nb = neighbor( seed, i );
1210 const int faces = subEntities( seed, 1 );
1211 for(
int face = 0; face<faces; ++ face )
1213 if( neighbor( nb, face ).equals(seed) )
1215 return indexInInside( nb, face );
1218 DUNE_THROW(InvalidStateException,
"inverse intersection not found");
1223 template <
class EntitySeed>
1225 outerNormal(
const EntitySeed& seed,
const int i )
const
1227 const int face = this->
template subEntitySeed<1>( seed, i ).index();
1228 const int normalIdx = face * GlobalCoordinate :: dimension ;
1229 GlobalCoordinate normal = copyToGlobalCoordinate( grid_.
face_normals + normalIdx );
1231 if( nb != seed.index() )
1238 template <
class EntitySeed>
1240 unitOuterNormal(
const EntitySeed& seed,
const int i )
const
1242 const int face = this->
template subEntitySeed<1>( seed, i ).index();
1243 if( seed.index() == grid_.
face_cells[ 2*face ] )
1245 return unitOuterNormals_[ face ];
1249 GlobalCoordinate normal = unitOuterNormals_[ face ];
1255 template <
class EntitySeed>
1256 GlobalCoordinate centroids(
const EntitySeed& seed )
const
1258 if( ! seed.isValid() )
1259 return GlobalCoordinate( 0 );
1261 const int index = GlobalCoordinate :: dimension * seed.index();
1262 const int codim = EntitySeed::codimension;
1263 assert( index >= 0 && index <
size( codim ) * GlobalCoordinate :: dimension );
1269 else if ( codim == 1 )
1273 else if( codim == dim )
1279 DUNE_THROW(InvalidStateException,
"codimension not implemented");
1280 return GlobalCoordinate( 0 );
1284 GlobalCoordinate copyToGlobalCoordinate(
const double* coords )
const
1286 GlobalCoordinate coordinate;
1287 for(
int i=0; i<GlobalCoordinate::dimension; ++i )
1289 coordinate[ i ] = coords[ i ];
1294 template <
class EntitySeed>
1295 double volumes(
const EntitySeed& seed )
const
1297 static const int codim = EntitySeed::codimension;
1298 if( codim == dim || ! seed.isValid() )
1304 assert( seed.isValid() );
1310 else if ( codim == 1 )
1316 DUNE_THROW(InvalidStateException,
"codimension not implemented");
1326 for(
int i=0; i<3; ++i )
1328 cartDims_[ i ] = grid_.
cartdims[ i ];
1332 const int numCells =
size( 0 );
1334 cellVertices_.resize( numCells );
1339 typedef std::array<int, 3> KeyType;
1340 std::map< const KeyType, const int > vertexFaceTags;
1341 const int vertexFacePattern [8][3] = {
1352 for(
int i=0; i<8; ++i )
1354 KeyType key; key.fill( 4 );
1355 for(
int j=0; j<dim; ++j )
1357 key[ j ] = vertexFacePattern[ i ][ j ];
1360 vertexFaceTags.insert( std::make_pair( key, i ) );
1363 for (
int c = 0; c < numCells; ++c)
1373 typedef std::map<int,int> vertexmap_t;
1374 typedef typename vertexmap_t :: iterator iterator;
1376 std::vector< vertexmap_t > cell_pts( dim*2 );
1385 const int node = grid_.
face_nodes[ nodepos ];
1386 iterator it = cell_pts[ faceTag ].find( node );
1387 if( it == cell_pts[ faceTag ].end() )
1389 cell_pts[ faceTag ].insert( std::make_pair( node, 1 ) );
1399 typedef std::map< int, std::set<int> > vertexlist_t;
1400 vertexlist_t vertexList;
1402 for(
int faceTag = 0; faceTag<dim*2; ++faceTag )
1404 for( iterator it = cell_pts[ faceTag ].begin(),
1405 end = cell_pts[ faceTag ].end(); it != end; ++it )
1409 if( (*it).second == 1 )
1411 vertexList[ (*it).first ].insert( faceTag );
1416 assert(
int(vertexList.size()) == ( dim == 2 ? 4 : 8) );
1418 cellVertices_[ c ].resize( vertexList.size() );
1419 for(
auto it = vertexList.begin(), end = vertexList.end(); it != end; ++it )
1421 assert( (*it).second.size() == dim );
1422 KeyType key; key.fill( 4 );
1424 std::copy( (*it).second.begin(), (*it).second.end(), key.begin() );
1425 auto vx = vertexFaceTags.find( key );
1426 assert( vx != vertexFaceTags.end() );
1427 if( vx != vertexFaceTags.end() )
1429 if( (*vx).second >=
int(cellVertices_[ c ].size()) )
1430 cellVertices_[ c ].resize( (*vx).second+1 );
1432 cellVertices_[ c ][ (*vx).second ] = (*it).first ;
1438 geomTypes_.resize(dim + 1);
1440 for (
int codim = 0; codim <= dim; ++codim)
1442 tmp = Dune::GeometryTypes::cube(dim - codim);
1443 geomTypes_[codim].push_back(tmp);
1449 int minVx = std::numeric_limits<int>::max();
1451 for (
int c = 0; c < numCells; ++c)
1453 std::set<int> cell_pts;
1459 cell_pts.insert(fnbeg, fnend);
1462 cellVertices_[ c ].resize( cell_pts.size() );
1463 std::copy(cell_pts.begin(), cell_pts.end(), cellVertices_[ c ].begin() );
1464 maxVx = std::max( maxVx,
int( cell_pts.size() ) );
1465 minVx = std::min( minVx,
int( cell_pts.size() ) );
1468 if( minVx == maxVx && maxVx == 4 )
1470 for (
int c = 0; c < numCells; ++c)
1472 assert( cellVertices_[ c ].
size() == 4 );
1473 GlobalCoordinate center( 0 );
1474 GlobalCoordinate p[ dim+1 ];
1475 for(
int i=0; i<dim+1; ++i )
1477 const int vertex = cellVertices_[ c ][ i ];
1479 for(
int d=0; d<dim; ++d )
1486 for(
int d=0; d<dim; ++d )
1491 Dune::GeometryType simplex;
1492 simplex = Dune::GeometryTypes::simplex(dim);
1494 typedef Dune::AffineGeometry< ctype, dim, dimworld> AffineGeometryType;
1495 AffineGeometryType geometry( simplex, p );
1503 for(
int face = 0 ; face < faces; ++face )
1506 const int b = grid_.
face_cells[ 2*face + 1 ];
1508 assert( a >=0 || b >=0 );
1513 GlobalCoordinate centerDiff( 0 );
1516 for(
int d=0; d<dimworld; ++d )
1523 for(
int d=0; d<dimworld; ++d )
1531 for(
int d=0; d<dimworld; ++d )
1538 for(
int d=0; d<dimworld; ++d )
1544 GlobalCoordinate normal( 0 );
1545 for(
int d=0; d<dimworld; ++d )
1550 if( centerDiff.two_norm() < 1e-10 )
1554 if( centerDiff * normal < 0 )
1562 bool allSimplex = true ;
1563 bool allCube = true ;
1565 for (
int c = 0; c < numCells; ++c)
1567 const int nVx = cellVertices_[ c ].size();
1579 geomTypes_.resize(dim + 1);
1581 for (
int codim = 0; codim <= dim; ++codim)
1585 tmp = Dune::GeometryTypes::simplex(dim - codim);
1586 geomTypes_[ codim ].push_back( tmp );
1590 tmp = Dune::GeometryTypes::cube(dim - codim);
1591 geomTypes_[ codim ].push_back( tmp );
1595 tmp = Dune::GeometryTypes::none(dim - codim);
1596 geomTypes_[ codim ].push_back( tmp );
1606 const int normalIdx = face * GlobalCoordinate :: dimension ;
1607 GlobalCoordinate normal = copyToGlobalCoordinate( grid_.
face_normals + normalIdx );
1608 normal /= normal.two_norm();
1609 unitOuterNormals_[ face ] = normal;
1611 if( isBoundaryFace( face ) )
1615 const int facePos = 2 * face ;
1620 grid_.
face_cells[ facePos ] = -nBndSegments_;
1622 else if ( grid_.
face_cells[ facePos+1 ] < 0 )
1624 grid_.
face_cells[ facePos+1 ] = -nBndSegments_;
1630 void print( std::ostream& out,
const UnstructuredGridType& grid )
const
1633 for(
int c=0; c<numCells; ++c )
1635 out <<
"cell " << c <<
" : faces = " << std::endl;
1641 out << f <<
" vx = " ;
1642 while( fnbeg != fnend )
1644 out << *fnbeg <<
" ";
1651 const auto& vx = cellVertices_[ c ];
1652 out <<
"cell " << c <<
" : vertices = ";
1653 for(
size_t i=0; i<vx.size(); ++i )
1654 out << vx[ i ] <<
" ";
1661 UnstructuredGridPtr gridPtr_;
1662 const UnstructuredGridType& grid_;
1664 CommunicationType comm_;
1665 std::array< int, 3 > cartDims_;
1666 std::vector< std::vector< GeometryType > > geomTypes_;
1667 std::vector< std::vector< int > > cellVertices_;
1669 std::vector< GlobalCoordinate > unitOuterNormals_;
1675 size_t nBndSegments_;
1687 template<
int dim,
int dimworld,
typename coord_t >
1688 template<
int codim >
1738 template< PartitionIteratorType pitype >
1742 ::template Partition< pitype >::LeafIterator
1745 ::template Partition< pitype >::LevelIterator
1772#include <opm/grid/polyhedralgrid/persistentcontainer.hh>
1773#include <opm/grid/polyhedralgrid/cartesianindexmapper.hh>
1774#include <opm/grid/polyhedralgrid/gridhelpers.hh>
Main OPM-Core grid data structure along with helper functions for construction, destruction and readi...
void destroy_grid(struct UnstructuredGrid *g)
Destroy and deallocate an UnstructuredGrid and all its data.
Definition: UnstructuredGrid.c:32
struct UnstructuredGrid * allocate_grid(size_t ndims, size_t ncells, size_t nfaces, size_t nfacenodes, size_t ncellfaces, size_t nnodes)
Allocate and initialise an UnstructuredGrid where pointers are set to location with correct size.
Definition: UnstructuredGrid.c:87
Routines to construct fully formed grid structures from a simple Cartesian (i.e., tensor product) des...
struct UnstructuredGrid * create_grid_cart2d(int nx, int ny, double dx, double dy)
Form geometrically Cartesian grid in two space dimensions with equally sized cells.
Definition: cart_grid.c:95
struct UnstructuredGrid * create_grid_hexa3d(int nx, int ny, int nz, double dx, double dy, double dz)
Form geometrically Cartesian grid in three space dimensions with equally sized cells.
Definition: cart_grid.c:59
Definition: entityseed.hh:16
Definition: entity.hh:152
Definition: geometry.hh:245
Definition: indexset.hh:24
Definition: intersectioniterator.hh:16
Definition: intersection.hh:22
Definition: iterator.hh:21
Definition: geometry.hh:264
identical grid wrapper
Definition: grid.hh:159
void postAdapt()
Definition: grid.hh:604
Traits::ctype ctype
type of vector coordinates (e.g., double)
Definition: grid.hh:308
int ghostSize(int, int codim) const
obtain size of ghost region for a grid level
Definition: grid.hh:644
bool loadBalance(CommDataHandleIF< DataHandle, Data > &)
rebalance the load each process has to handle
Definition: grid.hh:748
Partition< All_Partition >::LevelGridView LevelGridView
View types for All_Partition.
Definition: grid.hh:245
Traits::template Codim< EntitySeed::codimension >::EntityPointer entityPointer(const EntitySeed &seed) const
obtain EntityPointer from EntitySeed.
Definition: grid.hh:808
Traits::template Codim< EntitySeed::codimension >::Entity entity(const EntitySeed &seed) const
obtain EntityPointer from EntitySeed.
Definition: grid.hh:819
bool preAdapt()
Definition: grid.hh:582
int ghostSize(int codim) const
obtain size of ghost region for the leaf grid
Definition: grid.hh:624
bool loadBalance(DofManager &)
rebalance the load each process has to handle
Definition: grid.hh:768
Traits::GlobalIdSet GlobalIdSet
type of global id set
Definition: grid.hh:283
int maxLevel() const
obtain maximal grid level
Definition: grid.hh:414
Traits::LevelIndexSet LevelIndexSet
type of level index set
Definition: grid.hh:271
LevelGridView levelGridView(int) const
View for a grid level for All_Partition.
Definition: grid.hh:792
PolyhedralGrid(const std::vector< int > &n, const std::vector< double > &dx)
constructor
Definition: grid.hh:346
Partition< pitype >::LeafGridView leafGridView() const
View for the leaf grid.
Definition: grid.hh:784
Traits::HierarchicIterator HierarchicIterator
iterator over the grid hierarchy
Definition: grid.hh:225
int size(int, int codim) const
obtain number of entites on a level
Definition: grid.hh:427
Partition< pitype >::LevelGridView levelGridView(int) const
View for a grid level.
Definition: grid.hh:775
void update()
update grid caches
Definition: grid.hh:838
void switchToDistributedView()
Switch to the distributed view.
Definition: grid.hh:698
int overlapSize(int, int) const
obtain size of overlap region for a grid level
Definition: grid.hh:634
int size(int codim) const
obtain number of leaf entities
Definition: grid.hh:438
void scatterData(DataHandle &handle) const
Moves data from the global (all data on process) view to the distributed view.
Definition: grid.hh:883
Traits::LocalIdSet LocalIdSet
type of local id set
Definition: grid.hh:300
typename Traits::Communication Communication
communicator with all other processes having some part of the grid
Definition: grid.hh:311
const CommunicationType & comm() const
obtain CollectiveCommunication object
Definition: grid.hh:711
int overlapSize(int) const
obtain size of overlap region for the leaf grid
Definition: grid.hh:615
Traits::LevelIntersectionIterator LevelIntersectionIterator
iterator over intersections with other entities on the same level
Definition: grid.hh:229
GridFamily::Traits Traits
type of the grid traits
Definition: grid.hh:208
PolyhedralGrid(UnstructuredGridPtr &&gridPtr)
constructor
Definition: grid.hh:365
Traits::LeafIntersectionIterator LeafIntersectionIterator
iterator over intersections with other entities on the leaf level
Definition: grid.hh:227
bool adapt(DataHandle &)
Definition: grid.hh:598
int size(int, GeometryType type) const
obtain number of entites on a level
Definition: grid.hh:467
void communicate(DataHandle &, InterfaceType, CommunicationDirection) const
communicate information on leaf entities
Definition: grid.hh:684
size_t numBoundarySegments() const
obtain number of leaf entities
Definition: grid.hh:487
void switchToGlobalView()
Switch to the global view.
Definition: grid.hh:692
bool loadBalance()
rebalance the load each process has to handle
Definition: grid.hh:727
PolyhedralGrid(const UnstructuredGridType &grid)
constructor
Definition: grid.hh:384
LeafGridView leafGridView() const
View for the leaf grid for All_Partition.
Definition: grid.hh:799
int size(GeometryType type) const
returns the number of boundary segments within the macro grid
Definition: grid.hh:476
void communicate(DataHandle &, InterfaceType, CommunicationDirection, int) const
communicate information on a grid level
Definition: grid.hh:663
bool adapt()
Definition: grid.hh:588
Traits::LeafIndexSet LeafIndexSet
type of leaf index set
Definition: grid.hh:261
Transform a corner-point grid ZCORN field to account for MINPV processing.
Definition: MinpvProcessor.hpp:35
Routines to form a complete UnstructuredGrid from a corner-point specification.
struct UnstructuredGrid * create_grid_cornerpoint(const struct grdecl *in, double tol)
Construct grid representation from corner-point specification of a particular geological model.
Definition: cornerpoint_grid.c:164
void compute_geometry(struct UnstructuredGrid *g)
Compute derived geometric primitives in a grid.
Definition: cornerpoint_grid.c:137
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
Low-level corner-point processing routines and supporting data structures.
traits structure containing types for a codimension
Definition: grid.hh:1691
Partition< All_Partition >::LeafIterator LeafIterator
type of level iterator
Definition: grid.hh:1756
Traits::template Codim< codim >::Entity Entity
type of entity
Definition: grid.hh:1700
Traits::template Codim< codim >::EntityPointer EntityPointer
type of entity pointer
Definition: grid.hh:1706
Traits::template Codim< codim >::LocalGeometry LocalGeometry
type of local geometry
Definition: grid.hh:1731
Partition< All_Partition >::LevelIterator LevelIterator
type of leaf iterator
Definition: grid.hh:1765
Traits::template Codim< codim >::Geometry Geometry
type of world geometry
Definition: grid.hh:1721
Types for GridView.
Definition: grid.hh:239
Data structure for an unstructured grid, unstructured meaning that any cell may have an arbitrary num...
Definition: UnstructuredGrid.h:99
int * face_nodes
Contains for each face, the indices of its adjacent nodes.
Definition: UnstructuredGrid.h:121
int * face_nodepos
For a face f, face_nodepos[f] contains the starting index for f's nodes in the face_nodes array.
Definition: UnstructuredGrid.h:127
int number_of_cells
The number of cells in the grid.
Definition: UnstructuredGrid.h:109
double * face_areas
Exact or approximate face areas.
Definition: UnstructuredGrid.h:173
int * cell_faces
Contains for each cell, the indices of its adjacent faces.
Definition: UnstructuredGrid.h:146
int number_of_faces
The number of faces in the grid.
Definition: UnstructuredGrid.h:111
double * cell_centroids
Exact or approximate cell centroids, stored consecutively for each cell.
Definition: UnstructuredGrid.h:192
int * cell_facepos
For a cell c, cell_facepos[c] contains the starting index for c's faces in the cell_faces array.
Definition: UnstructuredGrid.h:152
int * cell_facetag
If non-null, this array contains a number for cell-face adjacency indicating the face's position with...
Definition: UnstructuredGrid.h:244
double * cell_volumes
Exact or approximate cell volumes.
Definition: UnstructuredGrid.h:197
int * face_cells
For a face f, face_cells[2*f] and face_cells[2*f + 1] contain the cell indices of the cells adjacent ...
Definition: UnstructuredGrid.h:138
double * face_centroids
Exact or approximate face centroids, stored consecutively for each face.
Definition: UnstructuredGrid.h:168
int number_of_nodes
The number of nodes in the grid.
Definition: UnstructuredGrid.h:113
int cartdims[3]
Contains the size of the logical cartesian structure (if any) of the grid.
Definition: UnstructuredGrid.h:227
int * global_cell
If non-null, this array contains the logical cartesian indices (in a lexicographic ordering) of each ...
Definition: UnstructuredGrid.h:214
double * face_normals
Exact or approximate face normals, stored consecutively for each face.
Definition: UnstructuredGrid.h:184
double * node_coordinates
Node coordinates, stored consecutively for each node.
Definition: UnstructuredGrid.h:160
Raw corner-point specification of a particular geological model.
Definition: preprocess.h:56
const double * coord
Pillar end-points.
Definition: preprocess.h:58
int dims[3]
Cartesian box dimensions.
Definition: preprocess.h:57
const double * zcorn
Corner-point depths.
Definition: preprocess.h:59
const int * actnum
Explicit "active" map.
Definition: preprocess.h:60