37#ifndef OPM_ENTITY_HEADER
38#define OPM_ENTITY_HEADER
40#include <dune/common/version.hh>
41#include <dune/geometry/type.hh>
42#include <dune/grid/common/gridenums.hh>
44#include "PartitionTypeIndicator.hpp"
45#include <opm/grid/cpgrid/DefaultGeometryPolicy.hpp>
49 const std::vector<std::array<int,3>>&,
50 const std::vector<std::array<int,3>>&,
51 const std::vector<std::array<int,3>>&,
52 const std::vector<std::string>&);
59template<
int,
int>
class Geometry;
60template<
int,PartitionIteratorType>
class Iterator;
61class IntersectionIterator;
62class HierarchicIterator;
64class LevelGlobalIdSet;
70class Entity :
public EntityRep<codim>
72 friend class LevelGlobalIdSet;
73 friend class GlobalIdSet;
74 friend class HierarchicIterator;
75 friend class CpGridData;
77 const std::vector<std::array<int,3>>&,
78 const std::vector<std::array<int,3>>&,
79 const std::vector<std::array<int,3>>&,
80 const std::vector<std::string>&);
85 enum { codimension = codim };
86 enum { dimension = 3 };
87 enum { mydimension = dimension - codimension };
88 enum { dimensionworld = 3 };
110 typedef double ctype;
129 :
EntityRep<codim>(entityrep), pgrid_(&grid)
135 :
EntityRep<codim>(index_arg, orientation_arg), pgrid_(&grid)
140 Entity(
int index_arg,
bool orientation_arg)
141 :
EntityRep<codim>(index_arg, orientation_arg), pgrid_()
198 return Dune::GeometryTypes::cube(3 - codim);
306 static constexpr std::array<int,8> in_father_reference_elem_corner_indices_ = {0,1,2,3,4,5,6,7};
314#include "Iterators.hpp"
315#include "Intersection.hpp"
323 static_assert(codim == 0,
"");
330 static_assert(codim == 0,
"");
337 static_assert(codim == 0,
"");
344 static_assert(codim == 0,
"");
367 return pgrid_->partition_type_indicator_->getPartitionType(*
this);
373 return pgrid_->partition_type_indicator_->getPartitionTypeWhenLgrs(*
this, lgrsOnDistributedGrid);
379#include <opm/grid/cpgrid/CpGridData.hpp>
390 else if ( codim == 0 ){
401 return pgrid_->geomVector<codim>()[*
this];
408 static_assert(codim == 0,
"");
411 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
413 }
else if (cc == 3) {
414 assert(i >= 0 && i < 8);
415 int corner_index = pgrid_->cell_to_point_[this->index()][i];
416 typename Codim<cc>::Entity se(*pgrid_, corner_index,
true);
420 OPM_THROW(std::runtime_error,
421 "No subentity exists of codimension " + std::to_string(cc));
431 Iter end = ileafend();
432 for (Iter it = ileafbegin(); it != end; ++it) {
433 if (it->boundary())
return true;
459 return pgrid_->leaf_to_level_cells_.empty()? pgrid_->level_ : pgrid_->leaf_to_level_cells_[this-> index()][0];
473 if (pgrid_ -> parent_to_children_cells_.empty()){
477 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1);
491 return (isLeaf() && hasFather());
498 const auto refinementMark = pgrid_ -> getMark(*
this);
499 return (refinementMark == 1);
505 if ((pgrid_ -> child_to_parent_cells_.empty()) || (pgrid_ -> child_to_parent_cells_[this->index()][0] == -1)){
516 if (this->hasFather()){
517 const int& coarser_level = pgrid_ -> child_to_parent_cells_[this->index()][0];
518 const int& parent_cell_index = pgrid_ -> child_to_parent_cells_[this->index()][1];
519 return Entity<0>( *((*(pgrid_ -> level_data_ptr_))[coarser_level].get()), parent_cell_index,
true);
522 OPM_THROW(std::logic_error,
"Entity has no father.");
529 if (!(this->hasFather())){
530 OPM_THROW(std::logic_error,
"Entity has no father.");
532 auto idx_in_parent_cell = pgrid_ -> cell_to_idxInParentCell_[this->index()];
533 if (idx_in_parent_cell !=-1) {
534 const auto& cells_per_dim = (*(pgrid_ -> level_data_ptr_))[this->level()] -> cells_per_dim_;
535 const auto& auxArr = pgrid_ -> getReferenceRefinedCorners(idx_in_parent_cell, cells_per_dim);
536 FieldVector<double, 3> corners_in_father_reference_elem_temp[8] =
537 { auxArr[0], auxArr[1], auxArr[2], auxArr[3], auxArr[4], auxArr[5], auxArr[6], auxArr[7]};
538 auto in_father_reference_elem_corners = std::make_shared<EntityVariable<cpgrid::Geometry<0, 3>, 3>>();
541 mutable_in_father_reference_elem_corners.assign(corners_in_father_reference_elem_temp,
542 corners_in_father_reference_elem_temp + 8);
544 Dune::FieldVector<double, 3> center_in_father_reference_elem = {0., 0.,0.};
545 for (
int corn = 0; corn < 8; ++corn) {
546 for (
int c = 0; c < 3; ++c)
548 center_in_father_reference_elem[c] += corners_in_father_reference_elem_temp[corn][c]/8.;
552 double volume_in_father_reference_elem = double(1)/(cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]);
555 in_father_reference_elem_corners, in_father_reference_elem_corner_indices_.data());
558 OPM_THROW(std::logic_error,
"Entity has no father.");
567 return this->father();
569 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
571 const int& levelElem = pgrid_->leaf_to_level_cells_[this->index()][0];
572 const int& levelElemIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
587 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
589 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
593 throw std::invalid_argument(
"The entity provided does not belong to the leaf grid view. ");
603 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
605 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
616 const auto& level_data = (*(pgrid_ -> level_data_ptr_))[level()].get();
618 return level_data -> global_cell_[getLevelElem().index()];
[ provides Dune::Grid ]
Definition CpGrid.hpp:201
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:147
Represents an entity of a given codim, with positive or negative orientation.
Definition PartitionTypeIndicator.hpp:47
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 PartitionTypeIndicator.hpp:46
bool isRegular() const
Refinement is not defined for CpGrid.
Definition Entity.hpp:179
bool mightVanish() const
Returns true, if entity might disappear during the next call to adapt().
Definition Entity.hpp:496
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition Entity.hpp:514
bool operator!=(const Entity &other) const
Inequality.
Definition Entity.hpp:152
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:385
LeafIntersectionIterator ileafend() const
End leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:342
int getLevelCartesianIdx() const
Get Cartesian Index in the level grid view where the Entity was born.
Definition Entity.hpp:614
HierarchicIterator hend(int) const
Iterator end over the children/beyond last child iterator.
Definition Entity.hpp:358
const Geometry & geometry() const
Return the geometry of the entity (does not depend on its orientation).
Definition Entity.hpp:399
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:503
bool isValid() const
isValid method for EntitySeed
Definition Entity.hpp:439
Entity< 0 > getOrigin() const
Returns (1) parent entity in the level-grid the parent cell was born, if the entity was born in any L...
Definition Entity.hpp:563
Entity< 0 > getEquivLevelElem() const
Get equivalent element on the level grid where the entity was born, if grid = leaf-grid-view....
Definition Entity.hpp:598
HierarchicIterator hbegin(int) const
Iterator begin over the children. [If requested, also over descendants more than one generation away....
Definition Entity.hpp:350
EntitySeed seed() const
Return an entity seed (light-weight entity).
Definition Entity.hpp:160
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition Entity.hpp:196
PartitionType partitionType() const
In serial run, the only partitionType() is InteriorEntity.
Definition Entity.hpp:365
Dune::cpgrid::Geometry< 3, 3 > geometryInFather() const
Return LocalGeometry representing the embedding of the entity into its father (when hasFather() is tr...
Definition Entity.hpp:527
PartitionType partitionTypeWhenLgrs(bool) const
For parallel run, the entity - for now - does not see the CpGrid therefore we pass a bool to make the...
Definition Entity.hpp:371
Entity()
Constructor taking a grid and an integer entity representation.
Definition Entity.hpp:122
LeafIntersectionIterator ileafbegin() const
Start leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:335
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition Entity.hpp:266
Entity(int index_arg, bool orientation_arg)
Constructor taking a entity index, and orientation.
Definition Entity.hpp:140
LevelIntersectionIterator ilevelend() const
End level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:328
Entity(const CpGridData &grid, EntityRep< codim > entityrep)
Constructor taking a grid and an entity representation.
Definition Entity.hpp:128
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition Entity.hpp:482
Codim< cc >::Entity subEntity(int i) const
Obtain subentity.
Entity< 0 > getLevelElem() const
To be invoked only for leaf-grid-view entities. Get equivalent element on the level grid the leaf-ent...
Definition Entity.hpp:582
LevelIntersectionIterator ilevelbegin() const
Start level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:321
int level() const
Return the level of the entity in the grid hierarchy. Level = 0 represents the coarsest grid.
Definition Entity.hpp:447
bool isLeaf() const
Check if the entity is in the leafview.
Definition Entity.hpp:471
bool operator==(const Entity &other) const
Equality.
Definition Entity.hpp:146
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition Entity.hpp:426
Entity(const CpGridData &grid, int index_arg, bool orientation_arg)
Constructor taking a grid, entity index, and orientation.
Definition Entity.hpp:134
This class encapsulates geometry for vertices, intersections, and cells.
Definition Geometry.hpp:75
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Intersection.hpp:278
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10