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;
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>&);
86 static constexpr int dimension = 3;
87 static constexpr int mydimension = dimension -
codimension;
88 static constexpr int 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_()
195 return Dune::GeometryTypes::cube(3 - codim);
294 int getIdxInParentCell()
const;
305#include "Iterators.hpp"
306#include "Intersection.hpp"
314 static_assert(codim == 0,
"");
321 static_assert(codim == 0,
"");
328 static_assert(codim == 0,
"");
335 static_assert(codim == 0,
"");
358 return pgrid_->partition_type_indicator_->getPartitionType(*
this);
365#include <opm/grid/cpgrid/CpGridData.hpp>
376 else if ( codim == 0 ){
387 return pgrid_->geomVector<codim>()[*
this];
394 static_assert(codim == 0,
"");
397 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
399 }
else if (cc == 3) {
400 assert(i >= 0 && i < 8);
401 int corner_index = pgrid_->cell_to_point_[this->index()][i];
402 typename Codim<cc>::Entity se(*pgrid_, corner_index,
true);
406 OPM_THROW(std::runtime_error,
407 "No subentity exists of codimension " + std::to_string(cc));
417 Iter end = ileafend();
418 for (Iter it = ileafbegin(); it != end; ++it) {
419 if (it->boundary())
return true;
445 return pgrid_->leaf_to_level_cells_.empty()? pgrid_->level_ : pgrid_->leaf_to_level_cells_[this-> index()][0];
459 if (pgrid_ -> parent_to_children_cells_.empty()){
463 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1);
477 return (isLeaf() && hasFather());
484 const auto refinementMark = pgrid_ -> getMark(*
this);
485 return (refinementMark == 1);
491 if ((pgrid_ -> child_to_parent_cells_.empty()) || (pgrid_ -> child_to_parent_cells_[this->index()][0] == -1)){
502 if (this->hasFather()){
503 const int& coarser_level = pgrid_ -> child_to_parent_cells_[this->index()][0];
504 const int& parent_cell_index = pgrid_ -> child_to_parent_cells_[this->index()][1];
505 return Entity<0>( *((*(pgrid_ -> level_data_ptr_))[coarser_level].get()), parent_cell_index,
true);
508 OPM_THROW(std::logic_error,
"Entity has no father.");
515 return pgrid_ -> cell_to_idxInParentCell_[this->index()];
522 if (!(this->hasFather())){
523 OPM_THROW(std::logic_error,
"Entity has no father.");
527 static constexpr std::array<int,8> in_father_reference_elem_corner_indices = {0,1,2,3,4,5,6,7};
531 auto idx_in_parent_cell = pgrid_ -> cell_to_idxInParentCell_[this->index()];
532 if (idx_in_parent_cell !=-1) {
533 const auto& cells_per_dim = (*(pgrid_ -> level_data_ptr_))[this->level()] -> cells_per_dim_;
534 const auto& auxArr = pgrid_ -> getReferenceRefinedCorners(idx_in_parent_cell, cells_per_dim);
535 FieldVector<double, 3> corners_in_father_reference_elem_temp[8] =
536 { auxArr[0], auxArr[1], auxArr[2], auxArr[3], auxArr[4], auxArr[5], auxArr[6], auxArr[7]};
537 auto in_father_reference_elem_corners = std::make_shared<EntityVariable<cpgrid::Geometry<0, 3>, 3>>();
540 mutable_in_father_reference_elem_corners.assign(corners_in_father_reference_elem_temp,
541 corners_in_father_reference_elem_temp + 8);
543 Dune::FieldVector<double, 3> center_in_father_reference_elem = {0., 0.,0.};
544 for (
int corn = 0; corn < 8; ++corn) {
545 for (
int c = 0; c < 3; ++c)
547 center_in_father_reference_elem[c] += corners_in_father_reference_elem_temp[corn][c]/8.;
551 double volume_in_father_reference_elem = double(1)/(cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]);
554 in_father_reference_elem_corners, in_father_reference_elem_corner_indices.data());
557 OPM_THROW(std::logic_error,
"Entity has no father.");
566 return this->father();
568 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
570 const int& levelElem = pgrid_->leaf_to_level_cells_[this->index()][0];
571 const int& levelElemIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
586 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
588 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
599 const auto& level_data = (*(pgrid_ -> level_data_ptr_))[level()].get();
600 return level_data -> global_cell_[getLevelElem().index()];
[ provides Dune::Grid ]
Definition CpGrid.hpp:198
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:138
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
bool mightVanish() const
Returns true, if entity might disappear during the next call to adapt().
Definition Entity.hpp:482
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition Entity.hpp:500
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:371
LeafIntersectionIterator ileafend() const
End leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:333
EntitySeed seed() const
Return an entity seed (light-weight entity).
Definition Entity.hpp:160
int getLevelCartesianIdx() const
Get Cartesian Index in the level grid view where the Entity was born.
Definition Entity.hpp:597
HierarchicIterator hend(int) const
Iterator end over the children/beyond last child iterator.
Definition Entity.hpp:349
const Geometry & geometry() const
Return the geometry of the entity (does not depend on its orientation).
Definition Entity.hpp:385
Entity(const CpGridData &grid, int index_arg, bool orientation_arg)
Constructor taking a grid, entity index, and orientation.
Definition Entity.hpp:134
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:489
bool isValid() const
isValid method for EntitySeed
Definition Entity.hpp:425
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:562
bool operator==(const Entity &other) const
Equality.
Definition Entity.hpp:146
bool isRegular() const
Refinement is not defined for CpGrid.
Definition Entity.hpp:179
static constexpr int codimension
Definition Entity.hpp:85
Entity()
Constructor taking a grid and an integer entity representation.
Definition Entity.hpp:122
HierarchicIterator hbegin(int) const
Iterator begin over the children. [If requested, also over descendants more than one generation away....
Definition Entity.hpp:341
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition Entity.hpp:263
PartitionType partitionType() const
In serial run, the only partitionType() is InteriorEntity.
Definition Entity.hpp:356
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:520
Entity(int index_arg, bool orientation_arg)
Constructor taking a entity index, and orientation.
Definition Entity.hpp:140
LeafIntersectionIterator ileafbegin() const
Start leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:326
Entity(const CpGridData &grid, EntityRep< codim > entityrep)
Constructor taking a grid and an entity representation.
Definition Entity.hpp:128
LevelIntersectionIterator ilevelend() const
End level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:319
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition Entity.hpp:468
Entity< 0 > getLevelElem() const
Get equivalent element on the level grid where the entity was born, if grid = leaf-grid-view....
Definition Entity.hpp:581
LevelIntersectionIterator ilevelbegin() const
Start level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:312
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition Entity.hpp:193
int level() const
Return the level of the entity in the grid hierarchy. Level = 0 represents the coarsest grid.
Definition Entity.hpp:433
bool isLeaf() const
Check if the entity is in the leafview.
Definition Entity.hpp:457
bool operator!=(const Entity &other) const
Inequality.
Definition Entity.hpp:152
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition Entity.hpp:412
Codim< cc >::Entity subEntity(int i) const
Obtain subentity.
This class encapsulates geometry for vertices, intersections, and cells.
Definition Geometry.hpp:75
The global id set for Dune.
Definition Indexsets.hpp:450
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Intersection.hpp:278
Definition Indexsets.hpp:350
The namespace Dune is the main namespace for all Dune code.
Definition CartesianIndexMapper.hpp:10