My Project
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//
10// $Date$
11//
12// $Revision$
13//
14//===========================================================================
15
16/*
17 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
18 Copyright 2009, 2010, 2022 Equinor ASA.
19
20 This file is part of The Open Porous Media project (OPM).
21
22 OPM is free software: you can redistribute it and/or modify
23 it under the terms of the GNU General Public License as published by
24 the Free Software Foundation, either version 3 of the License, or
25 (at your option) any later version.
26
27 OPM is distributed in the hope that it will be useful,
28 but WITHOUT ANY WARRANTY; without even the implied warranty of
29 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 GNU General Public License for more details.
31
32 You should have received a copy of the GNU General Public License
33 along with OPM. If not, see <http://www.gnu.org/licenses/>.
34*/
35
36#ifndef OPM_ENTITY_HEADER
37#define OPM_ENTITY_HEADER
38
39#include <dune/common/version.hh>
40#include <dune/geometry/type.hh>
41#include <dune/grid/common/gridenums.hh>
42
43#include "PartitionTypeIndicator.hpp"
44
45
46namespace Dune
47{
48 namespace cpgrid
49 {
50
51 template<int,int> class Geometry;
52 template<int,PartitionIteratorType> class Iterator;
53 class IntersectionIterator;
54 class HierarchicIterator;
55 class CpGridData;
56 class LevelGlobalIdSet;
57
58
62 template <int codim>
63 class Entity : public EntityRep<codim>
64 {
65 friend class LevelGlobalIdSet;
66 friend class GlobalIdSet;
67 friend class HierarchicIterator;
68 friend class CpGridData;
69
70 public:
73 enum { codimension = codim };
74 enum { dimension = 3 };
75 enum { mydimension = dimension - codimension };
76 enum { dimensionworld = 3 };
77
78 // the official DUNE names
79 typedef Entity EntitySeed;
80
84 template <int cd>
85 struct Codim
86 {
88 };
89
90
91 typedef cpgrid::Geometry<3-codim,3> Geometry;
92 typedef Geometry LocalGeometry;
93
97
98 typedef double ctype;
99
104// Entity(const CpGridData& grid, int entityrep)
105// : EntityRep<codim>(entityrep), pgrid_(&grid)
106// {
107// }
108
111 : EntityRep<codim>(), pgrid_( 0 )
112 {
113 }
114
116 Entity(const CpGridData& grid, EntityRep<codim> entityrep)
117 : EntityRep<codim>(entityrep), pgrid_(&grid)
118 {
119 }
120
122 Entity(const CpGridData& grid, int index_arg, bool orientation_arg)
123 : EntityRep<codim>(index_arg, orientation_arg), pgrid_(&grid)
124 {
125 }
126
128 bool operator==(const Entity& other) const
129 {
130 return EntityRep<codim>::operator==(other) && pgrid_ == other.pgrid_;
131 }
132
134 bool operator!=(const Entity& other) const
135 {
136 return !operator==(other);
137 }
138
143 {
144 return EntitySeed( impl() );
145 }
146
148 const Geometry& geometry() const;
149
151 int level() const;
152
158 bool isLeaf() const;
159
161 bool isRegular() const
162 {
163 return true;
164 }
165
168 PartitionType partitionType() const;
169
172 GeometryType type() const
173 {
174 return Dune::GeometryTypes::cube(3 - codim);
175 }
176
178 unsigned int subEntities ( const unsigned int cc ) const;
179
182 template <int cc>
183 typename Codim<cc>::Entity subEntity(int i) const;
184
187
190
193
196
197
200
203
205 bool isNew() const
206 {
207 return false;
208 }
209
211 bool mightVanish() const
212 {
213 return false;
214 }
215
221 bool hasFather() const;
222
227
234
240
241 // Mimic Dune entity wrapper
243 const Entity& impl() const
244 {
245 return *this;
246 }
247
248 Entity& impl()
249 {
250 return *this;
251 }
252
255 bool isValid () const;
256
257 protected:
258 const CpGridData* pgrid_;
259 };
260
261 } // namespace cpgrid
262} // namespace Dune
263
264// now we include the Iterators.hh We need to do this here because for hbegin/hend the compiler
265// needs to know the size of hierarchicIterator
266#include "Iterators.hpp"
267#include "Intersection.hpp"
268namespace Dune
269{
270 namespace cpgrid
271 {
272 template<int codim>
274 {
275 static_assert(codim == 0, "");
276 return LevelIntersectionIterator(*pgrid_, *this, false);
277 }
278
279 template<int codim>
281 {
282 static_assert(codim == 0, "");
283 return LevelIntersectionIterator(*pgrid_, *this, true);
284 }
285
286 template<int codim>
288 {
289 static_assert(codim == 0, "");
290 return LeafIntersectionIterator(*pgrid_, *this, false);
291 }
292
293 template<int codim>
295 {
296 static_assert(codim == 0, "");
297 return LeafIntersectionIterator(*pgrid_, *this, true);
298 }
299
300
301 template<int codim>
303 {
304 // Creates iterator with first child as target if there is one. Otherwise empty stack and target.
305 return HierarchicIterator(*this, maxLevel);
306 }
307
309 template<int codim>
311 {
312 // Creates iterator with empty stack and target.
313 return HierarchicIterator(maxLevel);
314 }
315
316 template <int codim>
317 PartitionType Entity<codim>::partitionType() const
318 {
319 return pgrid_->partition_type_indicator_->getPartitionType(*this);
320 }
321 } // namespace cpgrid
322} // namespace Dune
323
324
325#include <opm/grid/cpgrid/CpGridData.hpp>
326
327namespace Dune {
328namespace cpgrid {
329
330template<int codim>
331unsigned int Entity<codim>::subEntities ( const unsigned int cc ) const
332{
333 if (cc == codim) {
334 return 1;
335 }
336 else if ( codim == 0 ){ // Cell/element/Entity<0>
337 if ( cc == 3 ) { // Get number of corners of the element.
338 return 8;
339 }
340 }
341 return 0;
342}
343
344template <int codim>
346{
347 return pgrid_->geomVector<codim>()[*this];
348}
349
350template <int codim>
351template <int cc>
352typename Entity<codim>::template Codim<cc>::Entity Entity<codim>::subEntity(int i) const
353{
354 static_assert(codim == 0, "");
355 if (cc == 0) { // Cell/element/Entity<0>
356 assert(i == 0);
357 typename Codim<cc>::Entity se(*pgrid_, this->index(), this->orientation());
358 return se;
359 } else if (cc == 3) { // Corner/Entity<3>
360 assert(i >= 0 && i < 8);
361 int corner_index = pgrid_->cell_to_point_[this->index()][i];
362 typename Codim<cc>::Entity se(*pgrid_, corner_index, true);
363 return se;
364 }
365 else {
366 OPM_THROW(std::runtime_error,
367 "No subentity exists of codimension " + std::to_string(cc));
368 }
369}
370
371template <int codim>
373{
374 // Copied implementation from EntityDefaultImplementation,
375 // except for not checking LevelIntersectionIterators.
376 typedef LeafIntersectionIterator Iter;
377 Iter end = ileafend();
378 for (Iter it = ileafbegin(); it != end; ++it) {
379 if (it->boundary()) return true;
380 }
381 return false;
382}
383
384template <int codim>
386{
387 return pgrid_ ? EntityRep<codim>::index() < pgrid_->size(codim) : false;
388}
389
390
391// level() It simply returns the level of the entity in the grid hierarchy.
392template <int codim>
394{
395 // if distributed_data_ is not empty, level_data_ptr_ has size 1.
396 if ((*(pgrid_ -> level_data_ptr_)).size() == 1){
397 return 0; // when there is no refinenment, level_ is not automatically instantiated
398 }
399 if (pgrid_ == (*(pgrid_->level_data_ptr_)).back().get()) { // entity on the leafview -> get the level where it was born:
400 return pgrid_ -> leaf_to_level_cells_[this-> index()][0]; // leaf_to_level_cells_ leafIdx -> {level/LGR, cell index in LGR}
401 }
402 else {
403 return pgrid_-> level_;
404 }
405}
406
407
408// isLeaf()
409// - if distributed_data_ is empty: an element is a leaf <-> hbegin and hend return the same iterator. Then,
410// *cells from level 0 (coarse grid) that are not parents, are Leaf.
411// *cells from any LGR are Leaf, since they do not children (nested refinement not supported).
412// - if distrubuted_data_ is NOT empty: there may be children on a different process. Therefore,
413// isLeaf() returns true <-> the element is a leaf entity of the global refinement hierarchy. Equivalently,
414// it can be checked whether parent_to_children_cells_ is empty.
415template<int codim>
417{
418 if ((pgrid_ -> parent_to_children_cells_).empty()){ // LGR cells
419 return true;
420 }
421 else {
422 return (std::get<0>((pgrid_ -> parent_to_children_cells_)[this-> index()]) == -1); // Cells from GLOBAL, not involved in any LGR
423 }
424}
425
426template<int codim>
428{
429 if (pgrid_ -> child_to_parent_cells_.empty()){
430 return false;
431 }
432 else{
433 const auto& [level, parent] = pgrid_-> child_to_parent_cells_[this->index()]; // {-1,-1};
434 return (parent != -1);
435 }
436}
437
438template<int codim>
440{
441 if (!(this->hasFather())){
442 OPM_THROW(std::logic_error, "Entity has no father.");
443 }
444 const int& coarse_level = pgrid_ -> child_to_parent_cells_[this->index()][0];
445 const int& parent_index = pgrid_ -> child_to_parent_cells_[this->index()][1];
446 const auto& coarse_grid = (*(pgrid_ -> level_data_ptr_))[coarse_level].get();
447 return Entity<0>( *coarse_grid, parent_index, true);
448}
449
450template<int codim>
452{
453 if (!(this->hasFather())){
454 OPM_THROW(std::logic_error, "Entity has no father.");
455 }
456 else{
457#if 0
458 //
459 DefaultGeometryPolicy local_geometry;
460 std::array<int,8> localEntity_to_point;
461 std::array<int,8> allcorners_localEntity;
462 EntityVariableBase<cpgrid::Geometry<0,3>>& local_corners = local_geometry.geomVector<codim>();
463 local_corners.resize(8);
464 // Get IJK index of the entity.
465 std::array<int,3> eIJK;
466 pgrid_ -> getIJK(this->index(), eIJK);
467 // Get dimension of the grid.
468 const auto& grid_dim = pgrid_ -> logicalCartesianSize(); // {cells_per_dim[0]*patch_dim[0] ...[1], ...[2]}
469 // Get the local coordinates of the entity (in the reference unit cube).
470 const std::vector<FieldVector<double, 3>>& local_corners_temp = {
471 // corner '0'
472 { double(eIJK[0])/grid_dim[0], double(eIJK[1])/grid_dim[1], double(eIJK[2])/grid_dim[2] },
473 // corner '1'
474 { (double(eIJK[0])+1)/grid_dim[0], double(eIJK[1])/grid_dim[1], double(eIJK[2])/grid_dim[2] },
475 // corner '2'
476 { double(eIJK[0])/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1], double(eIJK[2])/grid_dim[2] },
477 // corner '3'
478 { (double(eIJK[0])+1)/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1], double(eIJK[2])/grid_dim[2] },
479 // corner '4'
480 { double(eIJK[0])/grid_dim[0], double(eIJK[1])/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] },
481 // corner '5'
482 { (double(eIJK[0])+1)/grid_dim[0], double(eIJK[1])/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] },
483 // corner '6'
484 { double(eIJK[0])/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] },
485 // corner '7'
486 { (double(eIJK[0])+1)/grid_dim[0], (double(eIJK[1])+1)/grid_dim[1], (double(eIJK[2])+1)/grid_dim[2] }};
487 // Compute the center of the 'local-entity'.
488 Dune::FieldVector<double, 3> local_center = {0., 0.,0.};
489 for (int corn = 0; corn < 8; ++corn) {
490 local_corners[corn] = local_corners_temp[corn];
491 local_center += local_corners[corn].center()/8.;
492 }
493 // Compute the volume of the 'local-entity'.
494 double local_volume = double(1)/(grid_dim[0]*grid_dim[1]*grid_dim[2]);
495 // Indices of 'all the corners', in this case, 0-7 (required to construct a Geometry<3,3> object).
496 allcorners_localEntity= {0,1,2,3,4,5,6,7};
497 // Create a pointer to the first element of "cellfiedPatch_to_point" (required to construct a Geometry<3,3> object).
498 const int* localEntity_indices_storage_ptr = &allcorners_localEntity[0];
499 // Construct (and return) the Geometry<3,3> of the 'cellified patch'.
500 return Dune::cpgrid::Geometry<3,3>(local_center, local_volume,
501 local_geometry.geomVector<codim>(), localEntity_indices_storage_ptr);
502#endif
503 OPM_THROW(std::logic_error, "geometryInFather() not implemented");
504 return {};
505 }
506}
507
508} // namespace cpgrid
509} // namespace Dune
510
511
512#endif // OPM_ENTITY_HEADER
Struct that hods all the data needed to represent a Cpgrid.
Definition: CpGridData.hpp:122
Definition: DefaultGeometryPolicy.hpp:53
const EntityVariable< cpgrid::Geometry< 3 - codim, 3 >, codim > & geomVector() const
Definition: DefaultGeometryPolicy.hpp:81
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:64
bool mightVanish() const
Returns true, if entity might disappear during the next call to adapt(). Dummy.
Definition: Entity.hpp:211
Entity< 0 > father() const
ONLY FOR CELLS (Entity<0>).
Definition: Entity.hpp:439
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:331
LeafIntersectionIterator ileafend() const
End leaf-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:294
EntitySeed seed() const
Return an entity seed (light-weight entity).
Definition: Entity.hpp:142
HierarchicIterator hend(int) const
Iterator end over the children/beyond last child iterator.
Definition: Entity.hpp:310
const Geometry & geometry() const
Return the geometry of the entity (does not depend on its orientation).
Definition: Entity.hpp:345
Entity(const CpGridData &grid, int index_arg, bool orientation_arg)
Constructor taking a grid, entity index, and orientation.
Definition: Entity.hpp:122
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:427
bool isValid() const
isValid method for EntitySeed
Definition: Entity.hpp:385
bool operator==(const Entity &other) const
Equality.
Definition: Entity.hpp:128
bool isRegular() const
Refinement is not defined for CpGrid.
Definition: Entity.hpp:161
Entity()
Constructor taking a grid and an integer entity representation.
Definition: Entity.hpp:110
HierarchicIterator hbegin(int) const
Iterator begin over the children. [If requested, also over descendants more than one generation away....
Definition: Entity.hpp:302
const Entity & impl() const
Access the actual implementation class behind Entity interface class.
Definition: Entity.hpp:243
PartitionType partitionType() const
For now, the grid is serial and the only partitionType() is InteriorEntity.
Definition: Entity.hpp:317
Dune::cpgrid::Geometry< 3, 3 > geometryInFather() const
Return LocalGeometry representing the embedding of the entity inti its father (when hasFather() is tr...
Definition: Entity.hpp:451
LeafIntersectionIterator ileafbegin() const
Start leaf-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:287
Entity(const CpGridData &grid, EntityRep< codim > entityrep)
Constructor taking a grid and an entity representation.
Definition: Entity.hpp:116
LevelIntersectionIterator ilevelend() const
End level-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:280
bool isNew() const
Returns true, if the entity has been created during the last call to adapt(). Dummy.
Definition: Entity.hpp:205
LevelIntersectionIterator ilevelbegin() const
Start level-iterator for the cell-cell intersections of this entity.
Definition: Entity.hpp:273
GeometryType type() const
Return marker object (GeometryType object) representing the reference element of the entity.
Definition: Entity.hpp:172
int level() const
Return the level of the entity in the grid hierarchy. Level = 0 represents the coarsest grid.
Definition: Entity.hpp:393
bool isLeaf() const
Check if the entity is in the leafview.
Definition: Entity.hpp:416
bool operator!=(const Entity &other) const
Inequality.
Definition: Entity.hpp:134
bool hasBoundaryIntersections() const
Returns true if any of my intersections are on the boundary.
Definition: Entity.hpp:372
Codim< cc >::Entity subEntity(int i) const
Obtain subentity.
This class encapsulates geometry for vertices, intersections, and cells.
Definition: Geometry.hpp:74
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:288
Definition: Indexsets.hpp:248
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
Definition: Entity.hpp:86