1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- 2 // vi: set et ts=4 sw=2 sts=2: 3 #ifndef DUNE_UGGRIDENTITY_HH 4 #define DUNE_UGGRIDENTITY_HH 5 6 /** \file 7 * \brief The UGGridEntity class and its specializations 8 */ 9 10 #include <memory> 11 12 #include <dune/geometry/referenceelements.hh> 13 14 #include <dune/grid/common/gridenums.hh> 15 #include <dune/grid/uggrid/uggridrenumberer.hh> 16 17 18 namespace Dune { 19 20 // Forward declarations 21 template<int dim> 22 class UGGrid; 23 template<int codim, class GridImp> 24 class UGGridEntitySeed; 25 template<int codim, PartitionIteratorType pitype, class GridImp> 26 class UGGridLevelIterator; 27 template<class GridImp> 28 class UGGridLevelIntersectionIterator; 29 template<class GridImp> 30 class UGGridLeafIntersectionIterator; 31 template<class GridImp> 32 class UGGridHierarchicIterator; 33 template <class GridType> 34 class GridFactory; 35 36 //********************************************************************** 37 // 38 // --UGGridEntity 39 // --Entity 40 // 41 /** \brief The implementation of entities in a UGGrid 42 \ingroup UGGrid 43 44 A Grid is a container of grid entities. An entity is parametrized by the codimension. 45 An entity of codimension c in dimension d is a d-c dimensional object. 46 47 */ 48 template<int codim, int dim, class GridImp> 49 class UGGridEntity 50 { 51 52 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 53 friend class UGGridLeafIterator; 54 55 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 56 friend class UGGridLevelIterator; 57 58 friend class UGGrid<dim>; 59 60 template <class GridImp_> 61 friend class UGGridLevelIndexSet; 62 63 template <class GridImp_> 64 friend class UGGridLeafIndexSet; 65 66 template <class GridImp_> 67 friend class UGGridIdSet; 68 69 friend class UGGridEntitySeed<codim, GridImp>; 70 71 typedef typename GridImp::ctype UGCtype; 72 73 typedef typename GridImp::Traits::template Codim<codim>::GeometryImpl GeometryImpl; 74 75 template <class GridType> 76 friend class GridFactory; 77 78 public: UGGridEntity()79 UGGridEntity() 80 : target_(nullptr) 81 , gridImp_(nullptr) 82 {} 83 84 /** \brief Copy constructor */ UGGridEntity(const UGGridEntity & other)85 UGGridEntity(const UGGridEntity& other) 86 : target_(other.target_) 87 , gridImp_(other.gridImp_) 88 { 89 if constexpr (codim==dim) 90 geo_ = other.geo_; 91 else 92 geo_ = std::make_unique<GeometryImpl>(*other.geo_); 93 } 94 UGGridEntity(typename UG_NS<dim>::template Entity<codim>::T * target,const GridImp * gridImp)95 UGGridEntity(typename UG_NS<dim>::template Entity<codim>::T* target, const GridImp* gridImp) 96 { 97 setToTarget(target,gridImp); 98 } 99 100 typedef typename GridImp::template Codim<codim>::Geometry Geometry; 101 102 /** \brief The type of UGGrid Entity seeds */ 103 typedef typename GridImp::Traits::template Codim<codim>::EntitySeed EntitySeed; 104 105 //! level of this entity level() const106 int level () const { 107 return UG_NS<dim>::myLevel(target_); 108 } 109 110 /** \brief Return the entity type identifier */ 111 GeometryType type() const; 112 113 /** \brief The partition type for parallel computing */ partitionType() const114 PartitionType partitionType () const 115 { 116 #ifndef ModelP 117 return InteriorEntity; 118 #else 119 if (UG_NS<dim>::Priority(target_) == UG_NS<dim>::PrioHGhost 120 || UG_NS<dim>::Priority(target_) == UG_NS<dim>::PrioVGhost 121 || UG_NS<dim>::Priority(target_) == UG_NS<dim>::PrioVHGhost) 122 return GhostEntity; 123 else if (hasBorderCopy()) 124 return BorderEntity; 125 else if (UG_NS<dim>::Priority(target_) == UG_NS<dim>::PrioMaster || UG_NS<dim>::Priority(target_) == UG_NS<dim>::PrioNone) 126 return InteriorEntity; 127 else 128 DUNE_THROW(GridError, "Unknown priority " << UG_NS<dim>::Priority(target_)); 129 #endif 130 } 131 132 /** \brief Get the partition type of each copy of a distributed entity 133 * 134 * This is a non-interface method, intended mainly for debugging and testing. 135 * 136 * \return Each pair in the return value contains a process number and the corresponding partition type 137 */ partitionTypes() const138 std::vector<std::pair<int,PartitionType> > partitionTypes () const 139 { 140 std::vector<std::pair<int,PartitionType> > result; 141 142 #ifndef ModelP 143 result.push_back(std::make_pair(0, InteriorEntity)); 144 #else 145 int *plist = UG_NS<dim>::DDD_InfoProcList(gridImp_->multigrid_->dddContext(), 146 UG_NS<dim>::ParHdr(target_)); 147 148 for (int i = 0; plist[i] >= 0; i += 2) 149 { 150 int rank = plist[i]; 151 auto priority = plist[i + 1]; 152 153 if (priority == UG_NS<dim>::PrioHGhost || priority == UG_NS<dim>::PrioVGhost || priority == UG_NS<dim>::PrioVHGhost) 154 result.push_back(std::make_pair(rank, GhostEntity)); 155 else 156 { 157 // The entity is not ghost. If it is (UG)PrioBorder somewhere, it is (Dune)Border. 158 // Otherwise it is (Dune)Interior. 159 bool hasBorderCopy = false; 160 for (int i = 0; plist[i] >= 0; i += 2) 161 if (plist[i + 1] == UG_NS<dim>::PrioBorder) 162 { 163 hasBorderCopy = true; 164 break; 165 } 166 167 result.push_back(std::make_pair(rank, (hasBorderCopy) ? BorderEntity : InteriorEntity)); 168 } 169 } 170 #endif 171 return result; 172 } 173 174 protected: 175 #ifdef ModelP 176 // \todo Unify with the following method hasBorderCopy() const177 bool hasBorderCopy() const 178 { 179 int *plist = UG_NS<dim>::DDD_InfoProcList( 180 gridImp_->multigrid_->dddContext(), 181 UG_NS<dim>::ParHdr(target_)); 182 for (int i = 0; plist[i] >= 0; i += 2) 183 if (plist[i + 1] == UG_NS<dim>::PrioBorder) 184 return true; 185 186 return false; 187 } 188 #endif 189 190 191 public: 192 193 //! geometry of this entity geometry() const194 Geometry geometry () const 195 { 196 if constexpr ((dim-codim)==0) 197 return Geometry( geo_ ); 198 else 199 return Geometry( *geo_ ); 200 } 201 202 /** \brief Get the seed corresponding to this entity */ seed() const203 EntitySeed seed () const { return EntitySeed( *this ); } 204 205 /** \brief Return the number of subEntities of codimension codim. 206 */ subEntities(unsigned int cd) const207 unsigned int subEntities (unsigned int cd) const 208 { 209 return referenceElement<UGCtype, dim-codim>(type()).size(cd-codim); 210 } 211 getTarget() const212 typename UG_NS<dim>::template Entity<codim>::T* getTarget() const 213 { 214 return target_; 215 } 216 217 //! equality equals(const UGGridEntity & other) const218 bool equals(const UGGridEntity& other) const 219 { 220 return getTarget() == other.getTarget(); 221 } 222 223 private: 224 /** \brief Set this entity to a particular UG entity */ setToTarget(typename UG_NS<dim>::template Entity<codim>::T * target,const GridImp * gridImp)225 void setToTarget(typename UG_NS<dim>::template Entity<codim>::T* target,const GridImp* gridImp) { 226 gridImp_ = gridImp; 227 target_ = target; 228 229 if constexpr ((dim-codim)==1) // Edge entity 230 { 231 // Obtain the corner coordinates from UG 232 UGCtype* cornerCoords[2*dim]; 233 UG_NS<dim>::Corner_Coordinates(target_, cornerCoords); 234 235 // convert to the type required by MultiLinearGeometry 236 std::vector<FieldVector<UGCtype, dim> > geometryCoords(2); 237 for (size_t i=0; i < 2; i++) 238 for (size_t j=0; j < dim; j++) 239 geometryCoords[i][j] = cornerCoords[i][j]; 240 241 geo_ = std::make_unique<GeometryImpl>(type(), geometryCoords); 242 } 243 else if constexpr ((dim-codim)==2) // Facet entity 244 { 245 // obtain the corner coordinates from UG 246 UGCtype* cornerCoords[4*dim]; 247 UG_NS<dim>::Corner_Coordinates(target_, cornerCoords); 248 249 // convert to the type required by MultiLinearGeometry 250 size_t numCorners = type().isTriangle() ? 3 : 4; 251 std::vector<FieldVector<UGCtype, dim> > geometryCoords(numCorners); 252 for(size_t i = 0; i < numCorners; i++) 253 for (size_t j = 0; j < dim; j++) 254 geometryCoords[UGGridRenumberer<dim-1>::verticesUGtoDUNE(i, type())][j] = cornerCoords[i][j]; 255 256 geo_ = std::make_unique<GeometryImpl>(type(), geometryCoords); 257 } 258 else 259 geo_.setToTarget(target); 260 } 261 262 // The geometric realization of the entity. It is a native UG type for vertices, 263 // and a dune-geometry MultiLinearGeometry for edges and facets. 264 // TODO: Maybe use MultiLinearGeometry for all these cases! 265 std::conditional_t<(dim-codim)==0, GeometryImpl, std::unique_ptr<GeometryImpl> > geo_; 266 267 /** \brief The UG object that represents this entity 268 * - 0d entities: node 269 * - 1d entities: edge 270 * - 2d entities in 3d grids: side vector 271 */ 272 typename UG_NS<dim>::template Entity<codim>::T* target_; 273 274 const GridImp* gridImp_; 275 }; 276 277 /** \brief Specialization for codim-0-entities. 278 * \ingroup UGGrid 279 * 280 * This class embodies the topological parts of elements of the grid. 281 * It has an extended interface compared to the general entity class. 282 * For example, Entities of codimension 0 allow to visit all neighbors. 283 * 284 * UGGrid only implements the cases dim==dimworld==2 and dim=dimworld==3. 285 */ 286 template<int dim, class GridImp> 287 class UGGridEntity<0,dim,GridImp> 288 { 289 friend class UGGrid<dim>; 290 friend class UGGridLeafIntersectionIterator <GridImp>; 291 friend class UGGridHierarchicIterator <GridImp>; 292 293 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 294 friend class UGGridLeafIterator; 295 296 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 297 friend class UGGridLevelIterator; 298 299 typedef typename GridImp::ctype UGCtype; 300 301 typedef typename GridImp::Traits::template Codim< 0 >::GeometryImpl GeometryImpl; 302 typedef typename GridImp::Traits::template Codim< 0 >::LocalGeometryImpl LocalGeometryImpl; 303 304 public: 305 typedef typename GridImp::template Codim<0>::Geometry Geometry; 306 typedef typename GridImp::template Codim<0>::LocalGeometry LocalGeometry; 307 308 //! The Iterator over neighbors on this level 309 typedef UGGridLevelIntersectionIterator<GridImp> LevelIntersectionIterator; 310 311 //! The Iterator over neighbors on the leaf level 312 typedef UGGridLeafIntersectionIterator<GridImp> LeafIntersectionIterator; 313 314 //! Iterator over descendants of the entity 315 typedef UGGridHierarchicIterator<GridImp> HierarchicIterator; 316 317 /** \brief The type of UGGrid Entity seeds */ 318 typedef typename GridImp::Traits::template Codim<0>::EntitySeed EntitySeed; 319 UGGridEntity()320 UGGridEntity() 321 : target_(nullptr) 322 , gridImp_(nullptr) 323 {} 324 UGGridEntity(typename UG_NS<dim>::Element * target,const GridImp * gridImp)325 UGGridEntity(typename UG_NS<dim>::Element* target, const GridImp* gridImp) 326 { 327 setToTarget(target,gridImp); 328 } 329 330 //! Level of this element level() const331 int level () const { 332 return UG_NS<dim>::myLevel(target_); 333 } 334 335 /** \brief Return the entity type identifier */ 336 GeometryType type() const; 337 338 /** \brief The partition type for parallel computing */ partitionType() const339 PartitionType partitionType () const { 340 #ifndef ModelP 341 return InteriorEntity; 342 #else 343 if (UG_NS<dim>::EPriority(target_) == UG_NS<dim>::PrioHGhost 344 || UG_NS<dim>::EPriority(target_) == UG_NS<dim>::PrioVGhost 345 || UG_NS<dim>::EPriority(target_) == UG_NS<dim>::PrioVHGhost) 346 return GhostEntity; 347 else 348 return InteriorEntity; 349 #endif 350 } 351 352 /** \brief Get the partition type of each copy of a distributed entity 353 * 354 * This is a non-interface method, intended mainly for debugging and testing. 355 * 356 * \return Each pair in the return value contains a process number and the corresponding partition type 357 */ partitionTypes() const358 std::vector<std::pair<int,PartitionType> > partitionTypes () const 359 { 360 std::vector<std::pair<int,PartitionType> > result; 361 362 #ifndef ModelP 363 result.push_back(std::make_pair(0, InteriorEntity)); 364 #else 365 int *plist = UG_NS<dim>::DDD_InfoProcList(gridImp_->multigrid_->dddContext(), 366 &target_->ge.ddd); 367 368 for (int i = 0; plist[i] >= 0; i += 2) 369 { 370 int rank = plist[i]; 371 auto priority = plist[i + 1]; 372 373 if (priority == UG_NS<dim>::PrioHGhost || priority == UG_NS<dim>::PrioVGhost || priority == UG_NS<dim>::PrioVHGhost) 374 result.push_back(std::make_pair(rank, GhostEntity)); 375 else 376 result.push_back(std::make_pair(rank, InteriorEntity)); 377 } 378 #endif 379 return result; 380 } 381 382 //! Geometry of this entity geometry() const383 Geometry geometry () const { return Geometry( geo_ ); } 384 385 /** \brief Get the seed corresponding to this entity */ seed() const386 EntitySeed seed () const { return EntitySeed( *this ); } 387 388 /** \brief Return the number of subEntities of codimension codim. 389 */ subEntities(unsigned int codim) const390 unsigned int subEntities (unsigned int codim) const 391 { 392 if (dim==3) { 393 394 switch (codim) { 395 case 0 : 396 return 1; 397 case 1 : 398 return UG_NS<dim>::Sides_Of_Elem(target_); 399 case 2 : 400 return UG_NS<dim>::Edges_Of_Elem(target_); 401 case 3 : 402 return UG_NS<dim>::Corners_Of_Elem(target_); 403 } 404 405 } else { 406 407 switch (codim) { 408 case 0 : 409 return 1; 410 case 1 : 411 return UG_NS<dim>::Edges_Of_Elem(target_); 412 case 2 : 413 return UG_NS<dim>::Corners_Of_Elem(target_); 414 default : 415 // do nothing for codim = 3 416 break; 417 } 418 419 } 420 DUNE_THROW(GridError, "You can't call UGGridEntity<0,dim>::subEntities " 421 << "with dim==" << dim << " and codim==" << codim << "!"); 422 } 423 424 425 /** \brief Provide access to sub entity i of given codimension. Entities 426 * are numbered 0 ... subEntities(cc)-1 427 */ 428 template<int cc> 429 typename GridImp::template Codim<cc>::Entity subEntity (int i) const; 430 431 /** \todo It would be faster to not use -1 as the end marker but 432 number of sides instead */ ileafbegin() const433 UGGridLeafIntersectionIterator<GridImp> ileafbegin () const { 434 return UGGridLeafIntersectionIterator<GridImp>(target_, (isLeaf()) ? 0 : UG_NS<dim>::Sides_Of_Elem(target_),gridImp_); 435 } 436 ilevelbegin() const437 UGGridLevelIntersectionIterator<GridImp> ilevelbegin () const { 438 return UGGridLevelIntersectionIterator<GridImp>(target_, 0, gridImp_); 439 } 440 441 //! Reference to one past the last leaf neighbor ileafend() const442 UGGridLeafIntersectionIterator<GridImp> ileafend () const { 443 return UGGridLeafIntersectionIterator<GridImp>(target_, UG_NS<dim>::Sides_Of_Elem(target_), gridImp_); 444 } 445 446 //! Reference to one past the last level neighbor ilevelend() const447 UGGridLevelIntersectionIterator<GridImp> ilevelend () const { 448 return UGGridLevelIntersectionIterator<GridImp>(target_, UG_NS<dim>::Sides_Of_Elem(target_),gridImp_); 449 } 450 451 //! returns true if Entity has NO children isLeaf() const452 bool isLeaf() const { 453 return UG_NS<dim>::isLeaf(target_); 454 } 455 456 //! returns true if element is of regular type isRegular() const457 bool isRegular() const { 458 return UG_NS<dim>::isRegular(target_); 459 } 460 461 /** \brief Returns true if the entity has intersections with the boundary 462 */ hasBoundaryIntersections() const463 bool hasBoundaryIntersections() const { 464 return UG_NS<dim>::isBoundaryElement(target_); 465 } 466 467 //! Inter-level access to father element on coarser grid. 468 //! Assumes that meshes are nested. father() const469 typename GridImp::template Codim<0>::Entity father () const { 470 return typename GridImp::template Codim<0>::Entity(UGGridEntity(UG_NS<dim>::EFather(target_),gridImp_)); 471 } 472 473 //! returns true if father entity exists hasFather() const474 bool hasFather () const 475 { 476 return UG_NS<dim>::EFather(target_) != nullptr; 477 } 478 479 /*! Location of this element relative to the reference element element of the father. 480 */ 481 LocalGeometry geometryInFather () const; 482 483 /*! Inter-level access to son elements on higher levels<=maxlevel. 484 This is provided for sparsely stored nested unstructured meshes. 485 Returns iterator to first son. 486 */ 487 UGGridHierarchicIterator<GridImp> hbegin (int maxlevel) const; 488 489 //! Returns iterator to one past the last son 490 UGGridHierarchicIterator<GridImp> hend (int maxlevel) const; 491 492 //*************************************************************** 493 // Interface for Adaptation 494 //*************************************************************** 495 496 //! returns true, if entity was refined during last adaptation cycle 497 bool isNew() const; 498 499 /*! \brief 500 returns true, if entity might be coarsened during next adaptation 501 cycle, which is true for entities that have been marked for 502 coarsening or for entities that are not regular (i.e. isRegular 503 returns false) */ 504 bool mightVanish() const; 505 506 //! 507 void setToTarget(typename UG_NS<dim>::Element* target, const GridImp* gridImp); 508 getTarget() const509 typename UG_NS<dim>::template Entity<0>::T* getTarget() const 510 { 511 return target_; 512 } 513 514 //! equality equals(const UGGridEntity & other) const515 bool equals(const UGGridEntity& other) const 516 { 517 return getTarget() == other.getTarget(); 518 } 519 520 //! the current geometry 521 GeometryImpl geo_; 522 523 /** \brief The corresponding UG-internal data structure */ 524 typename UG_NS<dim>::Element* target_; 525 526 /** \brief Pointer to the grid that we are part of. 527 * 528 * We need that only to hand it over to the intersections, 529 * which need it. 530 */ 531 const GridImp* gridImp_; 532 533 }; // end of UGGridEntity codim = 0 534 535 } // namespace Dune 536 537 #endif 538