My Project
Loading...
Searching...
No Matches
Entity.hpp
1//===========================================================================
2//
3// File: Entity.hpp
4//
5// Created: Fri May 29 20:26:48 2009
6//
7// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8// Brd Skaflestad <bard.skaflestad@sintef.no>
9// Antonella Ritorto <antonella.ritorto@opm-op.com>
10//
11// $Date$
12//
13// $Revision$
14//
15//===========================================================================
16
17/*
18 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
19 Copyright 2009, 2010, 2022 Equinor ASA.
20
21 This file is part of The Open Porous Media project (OPM).
22
23 OPM is free software: you can redistribute it and/or modify
24 it under the terms of the GNU General Public License as published by
25 the Free Software Foundation, either version 3 of the License, or
26 (at your option) any later version.
27
28 OPM is distributed in the hope that it will be useful,
29 but WITHOUT ANY WARRANTY; without even the implied warranty of
30 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
31 GNU General Public License for more details.
32
33 You should have received a copy of the GNU General Public License
34 along with OPM. If not, see <http://www.gnu.org/licenses/>.
35*/
36
37#ifndef OPM_ENTITY_HEADER
38#define OPM_ENTITY_HEADER
39
40#include <dune/common/version.hh>
41#include <dune/geometry/type.hh>
42#include <dune/grid/common/gridenums.hh>
43
44#include "PartitionTypeIndicator.hpp"
45#include <opm/grid/cpgrid/DefaultGeometryPolicy.hpp>
46
47// To be able to test local and global ids of vertices
48void refinePatch_and_check(Dune::CpGrid&,
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>&);
53
54namespace Dune
55{
56namespace cpgrid
57{
58
59template<int,int> class Geometry;
60template<int,PartitionIteratorType> class Iterator;
61class IntersectionIterator;
62class HierarchicIterator;
63class CpGridData;
64class LevelGlobalIdSet;
65
69template <int codim>
70class Entity : public EntityRep<codim>
71{
72 friend class LevelGlobalIdSet;
73 friend class GlobalIdSet;
74 friend class HierarchicIterator;
75 friend class CpGridData;
76 friend void ::refinePatch_and_check(Dune::CpGrid&,
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>&);
81
82public:
85 static constexpr int codimension = codim;
86 static constexpr int dimension = 3;
87 static constexpr int mydimension = dimension - codimension;
88 static constexpr int dimensionworld = 3;
89
90 // the official DUNE names
91 typedef Entity EntitySeed;
92
96 template <int cd>
97 struct Codim
98 {
100 };
101
102
103 typedef cpgrid::Geometry<3-codim,3> Geometry;
104 typedef Geometry LocalGeometry;
105
109
110 typedef double ctype;
111
116 // Entity(const CpGridData& grid, int entityrep)
117 // : EntityRep<codim>(entityrep), pgrid_(&grid)
118 // {
119 // }
120
123 : EntityRep<codim>(), pgrid_( 0 )
124 {
125 }
126
128 Entity(const CpGridData& grid, EntityRep<codim> entityrep)
129 : EntityRep<codim>(entityrep), pgrid_(&grid)
130 {
131 }
132
134 Entity(const CpGridData& grid, int index_arg, bool orientation_arg)
135 : EntityRep<codim>(index_arg, orientation_arg), pgrid_(&grid)
136 {
137 }
138
140 Entity(int index_arg, bool orientation_arg)
141 : EntityRep<codim>(index_arg, orientation_arg), pgrid_()
142 {
143 }
144
146 bool operator==(const Entity& other) const
147 {
148 return EntityRep<codim>::operator==(other) && pgrid_ == other.pgrid_;
149 }
150
152 bool operator!=(const Entity& other) const
153 {
154 return !operator==(other);
155 }
156
161 {
162 return EntitySeed( impl() );
163 }
164
166 const Geometry& geometry() const;
167
169 int level() const;
170
176 bool isLeaf() const;
177
179 bool isRegular() const
180 {
181 return true;
182 }
183
189 PartitionType partitionType() const;
190
193 GeometryType type() const
194 {
195 return Dune::GeometryTypes::cube(3 - codim);
196 }
197
199 unsigned int subEntities ( const unsigned int cc ) const;
200
203 template <int cc>
204 typename Codim<cc>::Entity subEntity(int i) const;
205
208
211
214
217
218
221
224
226 bool isNew() const;
227
234 bool mightVanish() const;
235
241 bool hasFather() const;
242
247
254
260
261 // Mimic Dune entity wrapper
263 const Entity& impl() const
264 {
265 return *this;
266 }
267
268 Entity& impl()
269 {
270 return *this;
271 }
272
275 bool isValid () const;
276
287
290
293
294 int getIdxInParentCell() const;
295
296protected:
297 const CpGridData* pgrid_;
298};
299
300} // namespace cpgrid
301} // namespace Dune
302
303// now we include the Iterators.hh We need to do this here because for hbegin/hend the compiler
304// needs to know the size of hierarchicIterator
305#include "Iterators.hpp"
306#include "Intersection.hpp"
307namespace Dune
308{
309namespace cpgrid
310{
311template<int codim>
313{
314 static_assert(codim == 0, "");
315 return LevelIntersectionIterator(*pgrid_, *this, false);
316}
317
318template<int codim>
320{
321 static_assert(codim == 0, "");
322 return LevelIntersectionIterator(*pgrid_, *this, true);
323}
324
325template<int codim>
327{
328 static_assert(codim == 0, "");
329 return LeafIntersectionIterator(*pgrid_, *this, false);
330}
331
332template<int codim>
334{
335 static_assert(codim == 0, "");
336 return LeafIntersectionIterator(*pgrid_, *this, true);
337}
338
339
340template<int codim>
342{
343 // Creates iterator with first child as target if there is one. Otherwise empty stack and target.
344 return HierarchicIterator(*this, maxLevel);
345}
346
348template<int codim>
350{
351 // Creates iterator with empty stack and target.
352 return HierarchicIterator(maxLevel);
353}
354
355template <int codim>
356PartitionType Entity<codim>::partitionType() const
357{
358 return pgrid_->partition_type_indicator_->getPartitionType(*this);
359}
360
361} // namespace cpgrid
362} // namespace Dune
363
364
365#include <opm/grid/cpgrid/CpGridData.hpp>
366
367namespace Dune {
368namespace cpgrid {
369
370template<int codim>
371unsigned int Entity<codim>::subEntities ( const unsigned int cc ) const
372{
373 if (cc == codim) {
374 return 1;
375 }
376 else if ( codim == 0 ){ // Cell/element/Entity<0>
377 if ( cc == 3 ) { // Get number of corners of the element.
378 return 8;
379 }
380 }
381 return 0;
382}
383
384template <int codim>
386{
387 return pgrid_->geomVector<codim>()[*this];
388}
389
390template <int codim>
391template <int cc>
392typename Entity<codim>::template Codim<cc>::Entity Entity<codim>::subEntity(int i) const
393{
394 static_assert(codim == 0, "");
395 if (cc == 0) { // Cell/element/Entity<0>
396 assert(i == 0);
397 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
398 return se;
399 } else if (cc == 3) { // Corner/Entity<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);
403 return se;
404 }
405 else {
406 OPM_THROW(std::runtime_error,
407 "No subentity exists of codimension " + std::to_string(cc));
408 }
409}
410
411template <int codim>
413{
414 // Copied implementation from EntityDefaultImplementation,
415 // except for not checking LevelIntersectionIterators.
416 typedef LeafIntersectionIterator Iter;
417 Iter end = ileafend();
418 for (Iter it = ileafbegin(); it != end; ++it) {
419 if (it->boundary()) return true;
420 }
421 return false;
422}
423
424template <int codim>
426{
427 return pgrid_ ? EntityRep<codim>::index() < pgrid_->size(codim) : false;
428}
429
430
431// level() It simply returns the level of the entity in the grid hierarchy.
432template <int codim>
434{
435 // Parallel and LGRs:
436 // If the grid has been distributed and - after that - LGRs have been added, then level_data_ptr_
437 // points at distributed_data_ which has size > 1.
438 //
439 // Serial and LGRs:
440 // If the grid is not distributed and there LGRs have been added, level_data_ptr_ points at data_ which
441 // has size > 1.
442 //
443 // - leaf_to_level_cells_ is non-empty only on the leaf grid view of a grid that has been refined.
444 // - level_ is set equal to zero when instantiating a pgrid, and rewriten when it corresponds to a refined level grid.
445 return pgrid_->leaf_to_level_cells_.empty()? pgrid_->level_ : pgrid_->leaf_to_level_cells_[this-> index()][0];
446}
447
448// isLeaf()
449// - if distributed_data_ is empty: an element is a leaf <-> hbegin and hend return the same iterator. Then,
450// *cells from level 0 (coarse grid) that are not parents, are Leaf.
451// *cells from any LGR are Leaf, since they do not children (nested refinement not supported).
452// - if distrubuted_data_ is NOT empty: there may be children on a different process. Therefore,
453// isLeaf() returns true <-> the element is a leaf entity of the global refinement hierarchy. Equivalently,
454// it can be checked whether parent_to_children_cells_ is empty.
455
456template<int codim>
458{
459 if (pgrid_ -> parent_to_children_cells_.empty()){ // LGR cells
460 return true;
461 }
462 else {
463 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1); // Cells from GLOBAL, not involved in any LGR
464 }
465}
466
467template<int codim>
469{
470 // WIP
471 // For an Entity "to be new" means that
472 // - it has a father
473 // - it has been created in the last call of adapt()
474 // To determine that, we track the level of the Entity and
475 // check if it was marked for refinement in the previos state
476 // of current_view_data_
477 return (isLeaf() && hasFather());
478}
479
480
481template<int codim>
483{
484 const auto refinementMark = pgrid_ -> getMark(*this);
485 return (refinementMark == 1);
486}
487
488template<int codim>
490{
491 if ((pgrid_ -> child_to_parent_cells_.empty()) || (pgrid_ -> child_to_parent_cells_[this->index()][0] == -1)){
492 return false;
493 }
494 else{
495 return true;
496 }
497}
498
499template<int codim>
501{
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);
506 }
507 else{
508 OPM_THROW(std::logic_error, "Entity has no father.");
509 }
510}
511
512template<int codim>
514{
515 return pgrid_ -> cell_to_idxInParentCell_[this->index()];
516}
517
518
519template<int codim>
521{
522 if (!(this->hasFather())){
523 OPM_THROW(std::logic_error, "Entity has no father.");
524 }
525
526 // Indices of corners in entity's geometry in father reference element.
527 static constexpr std::array<int,8> in_father_reference_elem_corner_indices = {0,1,2,3,4,5,6,7};
528 // 'static': The returned object Geometry<3,3> stores a pointer to in_father_reference_elem_corner_indices. Therefore,
529 // this variable is declared static to prolongate its lifetime beyond this function (static storage duration).
530
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>>();
538 EntityVariableBase<cpgrid::Geometry<0, 3>>& mutable_in_father_reference_elem_corners = *in_father_reference_elem_corners;
539 // Assign the corners. Make use of the fact that pointers behave like iterators.
540 mutable_in_father_reference_elem_corners.assign(corners_in_father_reference_elem_temp,
541 corners_in_father_reference_elem_temp + 8);
542 // Compute the center of the 'local-entity'.
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)
546 {
547 center_in_father_reference_elem[c] += corners_in_father_reference_elem_temp[corn][c]/8.;
548 }
549 }
550 // Compute the volume of the 'local-entity'.
551 double volume_in_father_reference_elem = double(1)/(cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]);
552 // Construct (and return) the Geometry<3,3> of 'child-cell in the reference element of its father (unit cube)'.
553 return Dune::cpgrid::Geometry<3,3>(center_in_father_reference_elem, volume_in_father_reference_elem,
554 in_father_reference_elem_corners, in_father_reference_elem_corner_indices.data());
555 }
556 else {
557 OPM_THROW(std::logic_error, "Entity has no father.");
558 }
559}
560
561template<int codim>
563{
564 if (hasFather())
565 {
566 return this->father();
567 }
568 if (!(pgrid_ -> leaf_to_level_cells_.empty())) // entity on the LeafGridView
569 { // leaf_to_level_cells_ [leaf idx] = { level where entity was born, cell idx in that level}
570 const int& levelElem = pgrid_->leaf_to_level_cells_[this->index()][0];
571 const int& levelElemIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
572 return Dune::cpgrid::Entity<0>( *((*(pgrid_ -> level_data_ptr_))[levelElem].get()), levelElemIdx, true);
573 }
574 else
575 {
576 return *this;
577 }
578}
579
580template<int codim>
582{
583 // Check that the element belongs to the leaf grid view
584 // This is needed to get the index of the element in the level it was born.
585 // leaf_to_level_cells_ [leaf idx] = {level where the entity was born, equivalent cell idx in that level}
586 if (!(pgrid_ -> leaf_to_level_cells_.empty())) // entity on the LeafGridView
587 {
588 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
589 return Dune::cpgrid::Entity<0>( *((*(pgrid_ -> level_data_ptr_))[this->level()].get()), entityLevelIdx, true);
590 }
591 else {
592 return *this;
593 }
594}
595
596template<int codim>
598{
599 const auto& level_data = (*(pgrid_ -> level_data_ptr_))[level()].get();
600 return level_data -> global_cell_[getLevelElem().index()];
601}
602
603} // namespace cpgrid
604} // namespace Dune
605
606
607#endif // OPM_ENTITY_HEADER
[ 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
Definition Entity.hpp:71
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
Definition Entity.hpp:98