36#ifndef OPM_ENTITY_HEADER
37#define OPM_ENTITY_HEADER
39#include <dune/common/version.hh>
40#include <dune/geometry/type.hh>
41#include <dune/grid/common/gridenums.hh>
43#include "PartitionTypeIndicator.hpp"
51 template<
int,
int>
class Geometry;
52 template<
int,PartitionIteratorType>
class Iterator;
53 class IntersectionIterator;
54 class HierarchicIterator;
56 class LevelGlobalIdSet;
73 enum { codimension = codim };
74 enum { dimension = 3 };
75 enum { mydimension = dimension - codimension };
76 enum { dimensionworld = 3 };
117 :
EntityRep<codim>(entityrep), pgrid_(&grid)
123 :
EntityRep<codim>(index_arg, orientation_arg), pgrid_(&grid)
174 return Dune::GeometryTypes::cube(3 - codim);
266#include "Iterators.hpp"
267#include "Intersection.hpp"
275 static_assert(codim == 0,
"");
282 static_assert(codim == 0,
"");
289 static_assert(codim == 0,
"");
296 static_assert(codim == 0,
"");
319 return pgrid_->partition_type_indicator_->getPartitionType(*
this);
325#include <opm/grid/cpgrid/CpGridData.hpp>
336 else if ( codim == 0 ){
347 return pgrid_->geomVector<codim>()[*
this];
354 static_assert(codim == 0,
"");
357 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
359 }
else if (cc == 3) {
360 assert(i >= 0 && i < 8);
361 int corner_index = pgrid_->cell_to_point_[this->index()][i];
362 typename Codim<cc>::Entity se(*pgrid_, corner_index,
true);
366 OPM_THROW(std::runtime_error,
367 "No subentity exists of codimension " + std::to_string(cc));
377 Iter end = ileafend();
378 for (Iter it = ileafbegin(); it != end; ++it) {
379 if (it->boundary())
return true;
396 if ((*(pgrid_ -> level_data_ptr_)).size() == 1){
399 if (pgrid_ == (*(pgrid_->level_data_ptr_)).back().get()) {
400 return pgrid_ -> leaf_to_level_cells_[this-> index()][0];
403 return pgrid_-> level_;
418 if ((pgrid_ -> parent_to_children_cells_).empty()){
422 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1);
429 if (pgrid_ -> child_to_parent_cells_.empty()){
433 const auto& [level, parent] = pgrid_-> child_to_parent_cells_[this->index()];
434 return (parent != -1);
441 if (!(this->hasFather())){
442 OPM_THROW(std::logic_error,
"Entity has no father.");
444 const int& coarse_level = pgrid_ -> child_to_parent_cells_[this->index()][0];
445 const int& parent_index = pgrid_ -> child_to_parent_cells_[this->index()][1];
446 const auto& coarse_grid = (*(pgrid_ -> level_data_ptr_))[coarse_level].get();
447 return Entity<0>( *coarse_grid, parent_index,
true);
453 if (!(this->hasFather())){
454 OPM_THROW(std::logic_error,
"Entity has no father.");
460 std::array<int,8> localEntity_to_point;
461 std::array<int,8> allcorners_localEntity;
463 local_corners.resize(8);
465 std::array<int,3> eIJK;
466 pgrid_ -> getIJK(this->index(), eIJK);
468 const auto& grid_dim = pgrid_ -> logicalCartesianSize();
470 const std::vector<FieldVector<double, 3>>& local_corners_temp = {
472 { double(eIJK[0])/grid_dim[0], double(eIJK[1])/grid_dim[1], double(eIJK[2])/grid_dim[2] },
474 { (double(eIJK[0])+1)/grid_dim[0],
double(eIJK[1])/grid_dim[1], double(eIJK[2])/grid_dim[2] },
476 { double(eIJK[0])/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1],
double(eIJK[2])/grid_dim[2] },
478 { (double(eIJK[0])+1)/grid_dim[0], (
double(eIJK[1])+1)/grid_dim[1],
double(eIJK[2])/grid_dim[2] },
480 { double(eIJK[0])/grid_dim[0], double(eIJK[1])/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] },
482 { (double(eIJK[0])+1)/grid_dim[0],
double(eIJK[1])/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] },
484 { double(eIJK[0])/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1], (
double(eIJK[2])+1)/grid_dim[2] },
486 { (double(eIJK[0])+1)/grid_dim[0], (
double(eIJK[1])+1)/grid_dim[1], (
double(eIJK[2])+1)/grid_dim[2] }};
488 Dune::FieldVector<double, 3> local_center = {0., 0.,0.};
489 for (
int corn = 0; corn < 8; ++corn) {
490 local_corners[corn] = local_corners_temp[corn];
491 local_center += local_corners[corn].center()/8.;
494 double local_volume = double(1)/(grid_dim[0]*grid_dim[1]*grid_dim[2]);
496 allcorners_localEntity= {0,1,2,3,4,5,6,7};
498 const int* localEntity_indices_storage_ptr = &allcorners_localEntity[0];
501 local_geometry.
geomVector<codim>(), localEntity_indices_storage_ptr);
503 OPM_THROW(std::logic_error,
"geometryInFather() not implemented");
Struct that hods all the data needed to represent a Cpgrid.
Definition: CpGridData.hpp:122
Definition: DefaultGeometryPolicy.hpp:53
const EntityVariable< cpgrid::Geometry< 3 - codim, 3 >, codim > & geomVector() const
Definition: DefaultGeometryPolicy.hpp:81
Represents an entity of a given codim, with positive or negative orientation.
Definition: EntityRep.hpp:99
bool operator==(const EntityRep &other) const
Equality operator.
Definition: EntityRep.hpp:179
int index() const
The (positive) index of an entity.
Definition: EntityRep.hpp:126
Base class for EntityVariable and SignedEntityVariable.
Definition: EntityRep.hpp:219
Definition: Entity.hpp:64
bool mightVanish() const
Returns true, if entity might disappear during the next call to adapt(). Dummy.
Definition: Entity.hpp:211
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition: Entity.hpp:439
unsigned int subEntities(const unsigned int cc) const
Return the number of all subentities of the entity of a given codimension cc.
Definition: Entity.hpp:331
LeafIntersectionIterator ileafend() const
End leaf-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:294
EntitySeed seed() const
Return an entity seed (light-weight entity).
Definition: Entity.hpp:142
HierarchicIterator hend(int) const
Iterator end over the children/beyond last child iterator.
Definition: Entity.hpp:310
const Geometry & geometry() const
Return the geometry of the entity (does not depend on its orientation).
Definition: Entity.hpp:345
Entity(const CpGridData &grid, int index_arg, bool orientation_arg)
Constructor taking a grid, entity index, and orientation.
Definition: Entity.hpp:122
bool hasFather() const
ONLY FOR CELLS (Entity<0>) Check if the entity comes from an LGR, i.e., it has been created via refin...
Definition: Entity.hpp:427
bool isValid() const
isValid method for EntitySeed
Definition: Entity.hpp:385
bool operator==(const Entity &other) const
Equality.
Definition: Entity.hpp:128
bool isRegular() const
Refinement is not defined for CpGrid.
Definition: Entity.hpp:161
Entity()
Constructor taking a grid and an integer entity representation.
Definition: Entity.hpp:110
HierarchicIterator hbegin(int) const
Iterator begin over the children. [If requested, also over descendants more than one generation away....
Definition: Entity.hpp:302
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition: Entity.hpp:243
PartitionType partitionType() const
For now, the grid is serial and the only partitionType() is InteriorEntity.
Definition: Entity.hpp:317
Dune::cpgrid::Geometry< 3, 3 > geometryInFather() const
Return LocalGeometry representing the embedding of the entity inti its father (when hasFather() is tr...
Definition: Entity.hpp:451
LeafIntersectionIterator ileafbegin() const
Start leaf-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:287
Entity(const CpGridData &grid, EntityRep< codim > entityrep)
Constructor taking a grid and an entity representation.
Definition: Entity.hpp:116
LevelIntersectionIterator ilevelend() const
End level-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:280
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition: Entity.hpp:205
LevelIntersectionIterator ilevelbegin() const
Start level-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:273
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition: Entity.hpp:172
int level() const
Return the level of the entity in the grid hierarchy. Level = 0 represents the coarsest grid.
Definition: Entity.hpp:393
bool isLeaf() const
Check if the entity is in the leafview.
Definition: Entity.hpp:416
bool operator!=(const Entity &other) const
Inequality.
Definition: Entity.hpp:134
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition: Entity.hpp:372
Codim< cc >::Entity subEntity(int i) const
Obtain subentity.
This class encapsulates geometry for vertices, intersections, and cells.
Definition: Geometry.hpp:74
The global id set for Dune.
Definition: Indexsets.hpp:304
Only needs to provide interface for doing nothing.
Definition: Iterators.hpp:108
Definition: Intersection.hpp:288
Definition: Indexsets.hpp:248
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
Definition: Entity.hpp:86