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 enum { codimension = codim };
86 enum { dimension = 3 };
87 enum { mydimension = dimension - codimension };
88 enum { 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
186 PartitionType partitionType() const;
187
192 PartitionType partitionTypeWhenLgrs(bool) const;
193
196 GeometryType type() const
197 {
198 return Dune::GeometryTypes::cube(3 - codim);
199 }
200
202 unsigned int subEntities ( const unsigned int cc ) const;
203
206 template <int cc>
207 typename Codim<cc>::Entity subEntity(int i) const;
208
211
214
217
220
221
224
227
229 bool isNew() const;
230
237 bool mightVanish() const;
238
244 bool hasFather() const;
245
250
257
263
264 // Mimic Dune entity wrapper
266 const Entity& impl() const
267 {
268 return *this;
269 }
270
271 Entity& impl()
272 {
273 return *this;
274 }
275
278 bool isValid () const;
279
290
293
296
299
300protected:
301 const CpGridData* pgrid_;
302private:
303 // static to not need any extra storage per Enitity. One object used for all instances
304 // constexpr to allow for in-class instantiation
306 static constexpr std::array<int,8> in_father_reference_elem_corner_indices_ = {0,1,2,3,4,5,6,7};
307};
308
309} // namespace cpgrid
310} // namespace Dune
311
312// now we include the Iterators.hh We need to do this here because for hbegin/hend the compiler
313// needs to know the size of hierarchicIterator
314#include "Iterators.hpp"
315#include "Intersection.hpp"
316namespace Dune
317{
318namespace cpgrid
319{
320template<int codim>
322{
323 static_assert(codim == 0, "");
324 return LevelIntersectionIterator(*pgrid_, *this, false);
325}
326
327template<int codim>
329{
330 static_assert(codim == 0, "");
331 return LevelIntersectionIterator(*pgrid_, *this, true);
332}
333
334template<int codim>
336{
337 static_assert(codim == 0, "");
338 return LeafIntersectionIterator(*pgrid_, *this, false);
339}
340
341template<int codim>
343{
344 static_assert(codim == 0, "");
345 return LeafIntersectionIterator(*pgrid_, *this, true);
346}
347
348
349template<int codim>
351{
352 // Creates iterator with first child as target if there is one. Otherwise empty stack and target.
353 return HierarchicIterator(*this, maxLevel);
354}
355
357template<int codim>
359{
360 // Creates iterator with empty stack and target.
361 return HierarchicIterator(maxLevel);
362}
363
364template <int codim>
365PartitionType Entity<codim>::partitionType() const
366{
367 return pgrid_->partition_type_indicator_->getPartitionType(*this);
368}
369
370template <int codim>
371PartitionType Entity<codim>::partitionTypeWhenLgrs(bool lgrsOnDistributedGrid) const
372{
373 return pgrid_->partition_type_indicator_->getPartitionTypeWhenLgrs(*this, lgrsOnDistributedGrid);
374}
375} // namespace cpgrid
376} // namespace Dune
377
378
379#include <opm/grid/cpgrid/CpGridData.hpp>
380
381namespace Dune {
382namespace cpgrid {
383
384template<int codim>
385unsigned int Entity<codim>::subEntities ( const unsigned int cc ) const
386{
387 if (cc == codim) {
388 return 1;
389 }
390 else if ( codim == 0 ){ // Cell/element/Entity<0>
391 if ( cc == 3 ) { // Get number of corners of the element.
392 return 8;
393 }
394 }
395 return 0;
396}
397
398template <int codim>
400{
401 return pgrid_->geomVector<codim>()[*this];
402}
403
404template <int codim>
405template <int cc>
406typename Entity<codim>::template Codim<cc>::Entity Entity<codim>::subEntity(int i) const
407{
408 static_assert(codim == 0, "");
409 if (cc == 0) { // Cell/element/Entity<0>
410 assert(i == 0);
411 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
412 return se;
413 } else if (cc == 3) { // Corner/Entity<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);
417 return se;
418 }
419 else {
420 OPM_THROW(std::runtime_error,
421 "No subentity exists of codimension " + std::to_string(cc));
422 }
423}
424
425template <int codim>
427{
428 // Copied implementation from EntityDefaultImplementation,
429 // except for not checking LevelIntersectionIterators.
430 typedef LeafIntersectionIterator Iter;
431 Iter end = ileafend();
432 for (Iter it = ileafbegin(); it != end; ++it) {
433 if (it->boundary()) return true;
434 }
435 return false;
436}
437
438template <int codim>
440{
441 return pgrid_ ? EntityRep<codim>::index() < pgrid_->size(codim) : false;
442}
443
444
445// level() It simply returns the level of the entity in the grid hierarchy.
446template <int codim>
448{
449 // Parallel and LGRs:
450 // If the grid has been distributed and - after that - LGRs have been added, then level_data_ptr_
451 // points at distributed_data_ which has size > 1.
452 //
453 // Serial and LGRs:
454 // If the grid is not distributed and there LGRs have been added, level_data_ptr_ points at data_ which
455 // has size > 1.
456 //
457 // - leaf_to_level_cells_ is non-empty only on the leaf grid view of a grid that has been refined.
458 // - level_ is set equal to zero when instantiating a pgrid, and rewriten when it corresponds to a refined level grid.
459 return pgrid_->leaf_to_level_cells_.empty()? pgrid_->level_ : pgrid_->leaf_to_level_cells_[this-> index()][0];
460}
461
462// isLeaf()
463// - if distributed_data_ is empty: an element is a leaf <-> hbegin and hend return the same iterator. Then,
464// *cells from level 0 (coarse grid) that are not parents, are Leaf.
465// *cells from any LGR are Leaf, since they do not children (nested refinement not supported).
466// - if distrubuted_data_ is NOT empty: there may be children on a different process. Therefore,
467// isLeaf() returns true <-> the element is a leaf entity of the global refinement hierarchy. Equivalently,
468// it can be checked whether parent_to_children_cells_ is empty.
469
470template<int codim>
472{
473 if (pgrid_ -> parent_to_children_cells_.empty()){ // LGR cells
474 return true;
475 }
476 else {
477 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1); // Cells from GLOBAL, not involved in any LGR
478 }
479}
480
481template<int codim>
483{
484 // WIP
485 // For an Entity "to be new" means that
486 // - it has a father
487 // - it has been created in the last call of adapt()
488 // To determine that, we track the level of the Entity and
489 // check if it was marked for refinement in the previos state
490 // of current_view_data_
491 return (isLeaf() && hasFather());
492}
493
494
495template<int codim>
497{
498 const auto refinementMark = pgrid_ -> getMark(*this);
499 return (refinementMark == 1);
500}
501
502template<int codim>
504{
505 if ((pgrid_ -> child_to_parent_cells_.empty()) || (pgrid_ -> child_to_parent_cells_[this->index()][0] == -1)){
506 return false;
507 }
508 else{
509 return true;
510 }
511}
512
513template<int codim>
515{
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);
520 }
521 else{
522 OPM_THROW(std::logic_error, "Entity has no father.");
523 }
524}
525
526template<int codim>
528{
529 if (!(this->hasFather())){
530 OPM_THROW(std::logic_error, "Entity has no father.");
531 }
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>>();
539 EntityVariableBase<cpgrid::Geometry<0, 3>>& mutable_in_father_reference_elem_corners = *in_father_reference_elem_corners;
540 // Assign the corners. Make use of the fact that pointers behave like iterators.
541 mutable_in_father_reference_elem_corners.assign(corners_in_father_reference_elem_temp,
542 corners_in_father_reference_elem_temp + 8);
543 // Compute the center of the 'local-entity'.
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)
547 {
548 center_in_father_reference_elem[c] += corners_in_father_reference_elem_temp[corn][c]/8.;
549 }
550 }
551 // Compute the volume of the 'local-entity'.
552 double volume_in_father_reference_elem = double(1)/(cells_per_dim[0]*cells_per_dim[1]*cells_per_dim[2]);
553 // Construct (and return) the Geometry<3,3> of 'child-cell in the reference element of its father (unit cube)'.
554 return Dune::cpgrid::Geometry<3,3>(center_in_father_reference_elem, volume_in_father_reference_elem,
555 in_father_reference_elem_corners, in_father_reference_elem_corner_indices_.data());
556 }
557 else {
558 OPM_THROW(std::logic_error, "Entity has no father.");
559 }
560}
561
562template<int codim>
564{
565 if (hasFather())
566 {
567 return this->father();
568 }
569 if (!(pgrid_ -> leaf_to_level_cells_.empty())) // entity on the LeafGridView
570 { // leaf_to_level_cells_ [leaf idx] = { level where entity was born, cell idx in that level}
571 const int& levelElem = pgrid_->leaf_to_level_cells_[this->index()][0];
572 const int& levelElemIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
573 return Dune::cpgrid::Entity<0>( *((*(pgrid_ -> level_data_ptr_))[levelElem].get()), levelElemIdx, true);
574 }
575 else
576 {
577 return *this;
578 }
579}
580
581template<int codim>
583{
584 // Check that the element belongs to the leaf grid view
585 // This is needed to get the index of the element in the level it was born.
586 // leaf_to_level_cells_ [leaf idx] = {level where the entity was born, equivalent cell idx in that level}
587 if (!(pgrid_ -> leaf_to_level_cells_.empty())) // entity on the LeafGridView
588 {
589 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
590 return Dune::cpgrid::Entity<0>( *((*(pgrid_ -> level_data_ptr_))[this->level()].get()), entityLevelIdx, true);
591 }
592 else {
593 throw std::invalid_argument("The entity provided does not belong to the leaf grid view. ");
594 }
595}
596
597template<int codim>
599{
600 // Check if the element belongs to the leaf grid view
601 // This is needed to get the index of the element in the level it was born.
602 // leaf_to_level_cells_ [leaf idx] = {level where the entity was born, equivalent cell idx in that level}
603 if (!(pgrid_ -> leaf_to_level_cells_.empty())) // entity on the LeafGridView
604 {
605 const int& entityLevelIdx = pgrid_->leaf_to_level_cells_[this->index()][1];
606 return Dune::cpgrid::Entity<0>( *((*(pgrid_ -> level_data_ptr_))[this->level()].get()), entityLevelIdx, true);
607 }
608 else {
609 return *this;
610 }
611}
612
613template<int codim>
615{
616 const auto& level_data = (*(pgrid_ -> level_data_ptr_))[level()].get();
617 // getLevelElem() throws when the entity does not belong to the leaf grid view.
618 return level_data -> global_cell_[getLevelElem().index()];
619}
620
621} // namespace cpgrid
622} // namespace Dune
623
624
625#endif // OPM_ENTITY_HEADER
[ 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 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:496
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition Entity.hpp:514
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
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: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
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: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
bool operator==(const Entity &other) const
Equality.
Definition Entity.hpp:146
bool isRegular() const
Refinement is not defined for CpGrid.
Definition Entity.hpp:179
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:350
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition Entity.hpp:266
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(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:335
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:328
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition Entity.hpp:482
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
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition Entity.hpp:196
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
Inequality.
Definition Entity.hpp:152
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition Entity.hpp:426
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:431
Only needs to provide interface for doing nothing.
Definition Iterators.hpp:108
Definition Intersection.hpp:278
Definition Indexsets.hpp:349
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10
Definition Entity.hpp:98