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>
53template<
int,
int>
class Geometry;
54template<
int,PartitionIteratorType>
class Iterator;
55class IntersectionIterator;
56class HierarchicIterator;
58class LevelGlobalIdSet;
74 enum { codimension = codim };
75 enum { dimension = 3 };
76 enum { mydimension = dimension - codimension };
77 enum { dimensionworld = 3 };
118 :
EntityRep<codim>(entityrep), pgrid_(&grid)
124 :
EntityRep<codim>(index_arg, orientation_arg), pgrid_(&grid)
129 Entity(
int index_arg,
bool orientation_arg)
130 :
EntityRep<codim>(index_arg, orientation_arg), pgrid_()
181 return Dune::GeometryTypes::cube(3 - codim);
281 static constexpr std::array<int,8> in_father_reference_elem_corner_indices_ = {0,1,2,3,4,5,6,7};
289#include "Iterators.hpp"
290#include "Intersection.hpp"
298 static_assert(codim == 0,
"");
305 static_assert(codim == 0,
"");
312 static_assert(codim == 0,
"");
319 static_assert(codim == 0,
"");
342 return pgrid_->partition_type_indicator_->getPartitionType(*
this);
348#include <opm/grid/cpgrid/CpGridData.hpp>
359 else if ( codim == 0 ){
370 return pgrid_->geomVector<codim>()[*
this];
377 static_assert(codim == 0,
"");
380 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
382 }
else if (cc == 3) {
383 assert(i >= 0 && i < 8);
384 int corner_index = pgrid_->cell_to_point_[this->index()][i];
385 typename Codim<cc>::Entity se(*pgrid_, corner_index,
true);
389 OPM_THROW(std::runtime_error,
390 "No subentity exists of codimension " + std::to_string(cc));
400 Iter end = ileafend();
401 for (Iter it = ileafbegin(); it != end; ++it) {
402 if (it->boundary())
return true;
419 if ((*(pgrid_ -> level_data_ptr_)).size() == 1){
422 if (pgrid_ == (*(pgrid_->level_data_ptr_)).back().get()) {
423 return pgrid_ -> leaf_to_level_cells_[this-> index()][0];
426 return pgrid_-> level_;
442 if ((pgrid_ -> parent_to_children_cells_).empty()){
446 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1);
453 if ((pgrid_ -> child_to_parent_cells_.empty()) || (pgrid_ -> child_to_parent_cells_[this->index()][0] == -1)){
464 if (this->hasFather()){
465 const int& coarse_level = pgrid_ -> child_to_parent_cells_[this->index()][0];
466 const int& parent_index = pgrid_ -> child_to_parent_cells_[this->index()][1];
467 const auto& coarse_grid = (*(pgrid_ -> level_data_ptr_))[coarse_level].get();
468 return Entity<0>( *coarse_grid, parent_index,
true);
471 OPM_THROW(std::logic_error,
"Entity has no father.");
478 if (!(this->hasFather())){
479 OPM_THROW(std::logic_error,
"Entity has no father.");
483 std::array<int,3> eIJK;
485 std::array<int,3> cells_per_dim;
487 std::array<int,3> child0_IJK;
488 const auto& level0_grid = (*(pgrid_->level_data_ptr_))[0];
489 const auto& child0_Idx = std::get<1>((*level0_grid).parent_to_children_cells_[this->father().index()])[0];
491 if (pgrid_ == (*(pgrid_->level_data_ptr_)).back().get())
493 const auto& lgr_grid = (*(pgrid_->level_data_ptr_))[
this -> level()];
494 cells_per_dim = (*lgr_grid).cells_per_dim_;
496 const auto& entity_lgrIdx = pgrid_ -> leaf_to_level_cells_[this->index()][1];
497 (*lgr_grid).getIJK(entity_lgrIdx, eIJK);
498 (*lgr_grid).getIJK(child0_Idx, child0_IJK);
502 pgrid_ -> getIJK(this->index(), eIJK);
503 cells_per_dim = pgrid_ -> cells_per_dim_;
504 pgrid_ -> getIJK(child0_Idx, child0_IJK);
509 FieldVector<double, 3> corners_in_father_reference_elem_temp[8] = {
511 { double(eIJK[0]-child0_IJK[0])/cells_per_dim[0], double(eIJK[1]-child0_IJK[1])/cells_per_dim[1],
512 double(eIJK[2]-child0_IJK[2])/cells_per_dim[2] },
514 { double(eIJK[0]-child0_IJK[0]+1)/cells_per_dim[0], double(eIJK[1]-child0_IJK[1])/cells_per_dim[1],
515 double(eIJK[2]-child0_IJK[2])/cells_per_dim[2] },
517 { double(eIJK[0]-child0_IJK[0])/cells_per_dim[0], double(eIJK[1]-child0_IJK[1]+1)/cells_per_dim[1],
518 double(eIJK[2]-child0_IJK[2])/cells_per_dim[2] },
520 { double(eIJK[0]-child0_IJK[0]+1)/cells_per_dim[0], double(eIJK[1]-child0_IJK[1]+1)/cells_per_dim[1],
521 double(eIJK[2]-child0_IJK[2])/cells_per_dim[2] },
523 { double(eIJK[0]-child0_IJK[0])/cells_per_dim[0], double(eIJK[1]-child0_IJK[1])/cells_per_dim[1],
524 double(eIJK[2]-child0_IJK[2]+1)/cells_per_dim[2] },
526 { double(eIJK[0]-child0_IJK[0]+1)/cells_per_dim[0], double(eIJK[1]-child0_IJK[1])/cells_per_dim[1],
527 double(eIJK[2]-child0_IJK[2]+1)/cells_per_dim[2] },
529 { double(eIJK[0]-child0_IJK[0])/cells_per_dim[0], double(eIJK[1]-child0_IJK[1]+1)/cells_per_dim[1],
530 double(eIJK[2]-child0_IJK[2]+1)/cells_per_dim[2] },
532 { double(eIJK[0]-child0_IJK[0]+1)/cells_per_dim[0], double(eIJK[1]-child0_IJK[1]+1)/cells_per_dim[1],
533 double(eIJK[2]-child0_IJK[2]+1)/cells_per_dim[2] }};
534 auto in_father_reference_elem_corners = std::make_shared<EntityVariable<cpgrid::Geometry<0, 3>, 3>>();
537 mutable_in_father_reference_elem_corners.assign(corners_in_father_reference_elem_temp,
538 corners_in_father_reference_elem_temp + 8);
540 Dune::FieldVector<double, 3> center_in_father_reference_elem = {0., 0.,0.};
541 for (
int corn = 0; corn < 8; ++corn) {
542 for (
int c = 0; c < 3; ++c)
544 center_in_father_reference_elem[c] += corners_in_father_reference_elem_temp[corn][c]/8.;
548 double volume_in_father_reference_elem = double(1)/(cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]);
551 in_father_reference_elem_corners, in_father_reference_elem_corner_indices_.data());
561 return this->father();
563 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
565 const int& entityIdxInLevel0 = pgrid_->leaf_to_level_cells_[this->index()][1];
566 const auto& coarse_grid = (*(pgrid_ -> level_data_ptr_))[0].get();
583 if (!(pgrid_ -> leaf_to_level_cells_.empty()))
585 const auto entityLevel =
this -> level();
586 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
587 const auto& level_grid = (*(pgrid_ -> level_data_ptr_))[entityLevel].get();
591 throw std::invalid_argument(
"The entity provided does not belong to the leaf grid view. ");
598 const auto entityLevel =
this -> level();
599 const auto level = (*(pgrid_ -> level_data_ptr_))[entityLevel].get();
600 const auto& elemInLevel = this->getLevelElem();
601 return level -> global_cell_[elemInLevel.index()];
Struct that hods all the data needed to represent a Cpgrid.
Definition CpGridData.hpp:131
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(). Dummy.
Definition Entity.hpp:218
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition Entity.hpp:462
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:354
LeafIntersectionIterator ileafend() const
End leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:317
EntitySeed seed() const
Return an entity seed (light-weight entity).
Definition Entity.hpp:149
int getLevelCartesianIdx() const
Get Cartesian Index in the level grid view where the Entity was born.
Definition Entity.hpp:596
HierarchicIterator hend(int) const
Iterator end over the children/beyond last child iterator.
Definition Entity.hpp:333
const Geometry & geometry() const
Return the geometry of the entity (does not depend on its orientation).
Definition Entity.hpp:368
Entity(const CpGridData &grid, int index_arg, bool orientation_arg)
Constructor taking a grid, entity index, and orientation.
Definition Entity.hpp:123
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:451
bool isValid() const
isValid method for EntitySeed
Definition Entity.hpp:408
Entity< 0 > getOrigin() const
getOrigin() Returns parent entity in level 0, if the entity was born in any LGR.
Definition Entity.hpp:557
bool operator==(const Entity &other) const
Equality.
Definition Entity.hpp:135
bool isRegular() const
Refinement is not defined for CpGrid.
Definition Entity.hpp:168
Entity()
Constructor taking a grid and an integer entity representation.
Definition Entity.hpp:111
HierarchicIterator hbegin(int) const
Iterator begin over the children. [If requested, also over descendants more than one generation away....
Definition Entity.hpp:325
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition Entity.hpp:250
PartitionType partitionType() const
For now, the grid is serial and the only partitionType() is InteriorEntity.
Definition Entity.hpp:340
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:476
Entity(int index_arg, bool orientation_arg)
Constructor taking a entity index, and orientation.
Definition Entity.hpp:129
LeafIntersectionIterator ileafbegin() const
Start leaf-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:310
Entity(const CpGridData &grid, EntityRep< codim > entityrep)
Constructor taking a grid and an entity representation.
Definition Entity.hpp:117
LevelIntersectionIterator ilevelend() const
End level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:303
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition Entity.hpp:212
Entity< 0 > getLevelElem() const
Get equivalent element on the level grid view for an element on the leaf grid view.
Definition Entity.hpp:576
LevelIntersectionIterator ilevelbegin() const
Start level-iterator for the cell-cell intersections of this entity.
Definition Entity.hpp:296
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition Entity.hpp:179
int level() const
Return the level of the entity in the grid hierarchy. Level = 0 represents the coarsest grid.
Definition Entity.hpp:416
bool isLeaf() const
Check if the entity is in the leafview.
Definition Entity.hpp:440
bool operator!=(const Entity &other) const
Inequality.
Definition Entity.hpp:141
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition Entity.hpp:395
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:304
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Intersection.hpp:278
Definition Indexsets.hpp:248
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10