1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- 2 // vi: set et ts=4 sw=2 sts=2: 3 4 #ifndef DUNE_UGGRID_HH 5 #define DUNE_UGGRID_HH 6 7 /** \file 8 * \brief The UGGrid class 9 */ 10 11 #include <memory> 12 13 #include <dune/common/classname.hh> 14 #include <dune/common/parallel/communication.hh> 15 #include <dune/common/exceptions.hh> 16 #include <dune/common/parallel/mpihelper.hh> 17 18 #include <dune/grid/common/boundarysegment.hh> 19 #include <dune/grid/common/capabilities.hh> 20 #include <dune/grid/common/grid.hh> 21 22 #if HAVE_UG || DOXYGEN 23 24 #ifdef ModelP 25 #include <dune/common/parallel/mpicommunication.hh> 26 #endif 27 28 /* [Before reading the following: the macros UG_DIM_2 and UG_DIM_3 where named 29 * _2 and _3, respectively, up until ug-3.12.0.] 30 * 31 * The following lines including the necessary UG headers are somewhat 32 tricky. Here's what's happening: 33 UG can support two- and three-dimensional grids. You choose by setting 34 either UG_DIM_2 or UG_DIM_3 while compiling. This changes all sorts of stuff, in 35 particular data structures in the headers. 36 UG was never supposed to provide 2d and 3d grids at the same time. 37 However, when compiling it as c++, the dimension-dependent parts are 38 wrapped up cleanly in the namespaces UG::D2 and UG::D3, respectively. That 39 way it is possible to link together the UG lib for 2d and the one for 3d. 40 But we also need the headers twice! Once with UG_DIM_2 set and once with UG_DIM_3! 41 So here we go:*/ 42 43 /* The following define tells the UG headers that we want access to a few 44 special fields, for example the extra index fields in the element data structures. 45 This define remains only for backwards compatibility with older version of UG. 46 All dune-uggrid versions since 2016-08-05 do not need this #define or the #undef 47 further below. */ 48 #define FOR_DUNE 49 50 // Set UG's space-dimension flag to 2d 51 #define UG_DIM_2 52 // And include all necessary UG headers 53 #include "uggrid/ugincludes.hh" 54 #undef DUNE_UGINCLUDES_HH 55 56 // Wrap a few large UG macros by functions before they get undef'ed away. 57 // Here: The 2d-version of the macros 58 #define UG_DIM 2 59 #include "uggrid/ugwrapper.hh" 60 #undef UG_DIM 61 62 // UG defines a whole load of preprocessor macros. ug_undefs.hh undefines 63 // them all, so we don't get name clashes. 64 #include "uggrid/ug_undefs.hh" 65 #undef UG_DIM_2 66 67 /* Now we're done with 2d, and we can do the whole thing over again for 3d */ 68 69 /* All macros set by UG have been unset. This includes the macros that ensure 70 single inclusion of headers. We can thus include them again. However, we 71 only want to include those headers again that contain dimension-dependent stuff. 72 Therefore, we set a few single-inclusion defines manually before including 73 ugincludes.hh again. 74 */ 75 #define UGTYPES_H 76 #define __HEAPS__ 77 #define __UGENV__ 78 #define __DEVICESH__ 79 #ifdef ModelP 80 #define __PPIF__ 81 #endif 82 83 #define UG_DIM_3 84 #include "uggrid/ugincludes.hh" 85 #undef DUNE_UGINCLUDES_HH 86 87 // Wrap a few large UG macros by functions before they get undef'ed away. 88 // This time it's the 3d-versions. 89 #define UG_DIM 3 90 #include "uggrid/ugwrapper.hh" 91 #undef UG_DIM 92 93 // undef all macros defined by UG 94 #include "uggrid/ug_undefs.hh" 95 96 #undef UG_DIM_3 97 #undef FOR_DUNE 98 99 // The components of the UGGrid interface 100 #include "uggrid/uggridgeometry.hh" 101 #include "uggrid/uggridlocalgeometry.hh" 102 #include "uggrid/uggridentity.hh" 103 #include "uggrid/uggridentityseed.hh" 104 #include "uggrid/uggridintersections.hh" 105 #include "uggrid/uggridintersectioniterators.hh" 106 #include "uggrid/uggridleveliterator.hh" 107 #include "uggrid/uggridleafiterator.hh" 108 #include "uggrid/uggridhieriterator.hh" 109 #include "uggrid/uggridindexsets.hh" 110 #include <dune/grid/uggrid/uggridviews.hh> 111 #ifdef ModelP 112 #include "uggrid/ugmessagebuffer.hh" 113 #include "uggrid/uglbgatherscatter.hh" 114 #endif 115 116 // Not needed here, but included for user convenience 117 #include "uggrid/uggridfactory.hh" 118 119 #ifdef ModelP 120 template <class DataHandle, int GridDim, int codim> 121 const Dune::UGGrid<GridDim>* Dune::UGMessageBuffer<DataHandle, GridDim, codim>::grid_; 122 123 template <class DataHandle, int GridDim, int codim> 124 DataHandle *Dune::UGMessageBuffer<DataHandle,GridDim,codim>::duneDataHandle_ = nullptr; 125 126 template <class DataHandle, int GridDim, int codim> 127 int Dune::UGMessageBuffer<DataHandle,GridDim,codim>::level = -1; 128 #endif // ModelP 129 130 namespace Dune { 131 132 #ifdef ModelP 133 using UGCollectiveCommunication = CollectiveCommunication<MPI_Comm>; 134 #else 135 using UGCollectiveCommunication = CollectiveCommunication<No_Comm>; 136 #endif 137 138 template<int dim> 139 struct UGGridFamily 140 { 141 typedef GridTraits<dim,dim,Dune::UGGrid<dim>, 142 UGGridGeometry, 143 UGGridEntity, 144 UGGridLevelIterator, 145 UGGridLeafIntersection, 146 UGGridLevelIntersection, 147 UGGridLeafIntersectionIterator, 148 UGGridLevelIntersectionIterator, 149 UGGridHierarchicIterator, 150 UGGridLeafIterator, 151 UGGridLevelIndexSet< const UGGrid<dim> >, 152 UGGridLeafIndexSet< const UGGrid<dim> >, 153 UGGridIdSet< const UGGrid<dim> >, 154 typename UG_NS<dim>::UG_ID_TYPE, 155 UGGridIdSet< const UGGrid<dim> >, 156 typename UG_NS<dim>::UG_ID_TYPE, 157 UGCollectiveCommunication, 158 UGGridLevelGridViewTraits, 159 UGGridLeafGridViewTraits, 160 UGGridEntitySeed, 161 UGGridLocalGeometry> 162 Traits; 163 }; 164 165 166 //********************************************************************** 167 // 168 // --UGGrid 169 // 170 //********************************************************************** 171 172 /** 173 \brief Front-end for the grid manager of the finite element toolbox UG3 174 175 \ingroup GridImplementations 176 177 This is the implementation of the grid interface using the UG3 grid management 178 system. It is best described in this <a href="https://doi.org/10.1007/s007910050003">paper</a>. 179 To our knowledge, the original code is not available anymore, 180 but the relevant parts have been forked into the %Dune module 181 dune-uggrid, available from <a href="https://www.dune-project.org/modules/dune-uggrid"></a>. 182 183 UGGrid provides conforming grids 184 in two and three space dimensions. The grids can be mixed, i.e. 185 2d grids can contain triangles and quadrilaterals and 3d grids can 186 contain tetrahedra and hexahedra and also pyramids and prisms. 187 The grid refinement rules are very flexible. Local adaptive red/green 188 refinement is the default, but a special method in the UGGrid class 189 allows you to directly access a number of anisotropic refinements rules. 190 Last but not least, the UG grid manager is completely parallelized, 191 and you can use boundaries parametrized by either analytical expressions 192 or high-resolution piecewise linear surfaces. 193 194 In your %Dune application, you can instantiate objects of the 195 type UGGrid<2> or UGGrid<3>. You can have more than one, if 196 you choose. It is even possible to have 2d and 3d grids at the same 197 time, even though the original UG system never intended to support 198 this! 199 200 See the documentation for the factory class GridFactory<UGGrid<dimworld> > 201 to learn how to create UGGrid objects. 202 */ 203 template <int dim> 204 class UGGrid : public GridDefaultImplementation <dim, dim, double, UGGridFamily<dim> > 205 { 206 typedef GridDefaultImplementation <dim, dim, double, UGGridFamily<dim> > Base; 207 208 friend class UGGridGeometry<0,dim,const UGGrid<dim> >; 209 friend class UGGridGeometry<dim,dim,const UGGrid<dim> >; 210 friend class UGGridGeometry<1,2,const UGGrid<dim> >; 211 friend class UGGridGeometry<2,3,const UGGrid<dim> >; 212 213 friend class UGGridEntity <0,dim,const UGGrid<dim> >; 214 friend class UGGridEntity <1,dim,const UGGrid<dim> >; 215 friend class UGGridEntity <2,dim,const UGGrid<dim> >; 216 friend class UGGridEntity <dim,dim,const UGGrid<dim> >; 217 friend class UGGridHierarchicIterator<const UGGrid<dim> >; 218 friend class UGGridLeafIntersection<const UGGrid<dim> >; 219 friend class UGGridLevelIntersection<const UGGrid<dim> >; 220 friend class UGGridLeafIntersectionIterator<const UGGrid<dim> >; 221 friend class UGGridLevelIntersectionIterator<const UGGrid<dim> >; 222 223 friend class UGGridLevelIndexSet<const UGGrid<dim> >; 224 friend class UGGridLeafIndexSet<const UGGrid<dim> >; 225 friend class UGGridIdSet<const UGGrid<dim> >; 226 template <class GridImp_> 227 friend class UGGridLeafGridView; 228 template <class GridImp_> 229 friend class UGGridLevelGridView; 230 231 friend class GridFactory<UGGrid<dim> >; 232 233 #ifdef ModelP 234 friend class UGLBGatherScatter; 235 #endif 236 237 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 238 friend class UGGridLeafIterator; 239 template <int codim_, PartitionIteratorType PiType_, class GridImp_> 240 friend class UGGridLevelIterator; 241 242 /** \brief UGGrid is only implemented for 2 and 3 dimension */ 243 static_assert(dim==2 || dim==3, "Use UGGrid only for 2d and 3d!"); 244 245 // The different instantiations are mutual friends so they can access 246 // each others numOfUGGrids field 247 friend class UGGrid<2>; 248 friend class UGGrid<3>; 249 250 //********************************************************** 251 // The Interface Methods 252 //********************************************************** 253 public: 254 //! type of the used GridFamily for this grid 255 typedef UGGridFamily<dim> GridFamily; 256 257 // the Traits 258 typedef typename UGGridFamily<dim>::Traits Traits; 259 260 //! The type used to store coordinates 261 typedef UG::DOUBLE ctype; 262 263 //! The type used for process ranks 264 typedef unsigned int Rank; 265 266 /** \brief Default constructor 267 */ 268 UGGrid(UGCollectiveCommunication comm = {}); 269 270 //! Destructor 271 ~UGGrid() noexcept(false); 272 273 //! Return maximum level defined in this grid. Levels are numbered 274 //! 0 ... maxlevel with 0 the coarsest level. 275 int maxLevel() const; 276 277 /** \brief Create an Entity from an EntitySeed */ 278 template <typename Seed> 279 typename Traits::template Codim<Seed::codimension>::Entity entity(const Seed & seed) const280 entity(const Seed& seed) const 281 { 282 const int codim = Seed::codimension; 283 return typename Traits::template Codim<codim>::Entity(UGGridEntity<codim,dim,const UGGrid<dim> >(seed.impl().target(),this)); 284 } 285 286 /** \brief Number of grid entities per level and codim 287 */ 288 int size (int level, int codim) const; 289 290 //! number of leaf entities per codim in this process size(int codim) const291 int size (int codim) const 292 { 293 return leafIndexSet().size(codim); 294 } 295 296 //! number of entities per level and geometry type in this process size(int level,GeometryType type) const297 int size (int level, GeometryType type) const 298 { 299 return this->levelIndexSet(level).size(type); 300 } 301 302 //! number of leaf entities per geometry type in this process size(GeometryType type) const303 int size (GeometryType type) const 304 { 305 return this->leafIndexSet().size(type); 306 } 307 308 /** \brief Return the number of boundary segments */ numBoundarySegments() const309 size_t numBoundarySegments() const { 310 // The number is stored as a member of UGGrid upon grid creation. 311 // The corresponding data structure is not exported by UG. (It is in ug/dom/std/std_internal.h) 312 return numBoundarySegments_; 313 } 314 315 /** \brief Access to the GlobalIdSet */ globalIdSet() const316 const typename Traits::GlobalIdSet& globalIdSet() const 317 { 318 return idSet_; 319 } 320 321 /** \brief Access to the LocalIdSet */ localIdSet() const322 const typename Traits::LocalIdSet& localIdSet() const 323 { 324 return idSet_; 325 } 326 327 /** \brief Access to the LevelIndexSets */ levelIndexSet(int level) const328 const typename Traits::LevelIndexSet& levelIndexSet(int level) const 329 { 330 if (level<0 || level>maxLevel()) 331 DUNE_THROW(GridError, "levelIndexSet of nonexisting level " << level << " requested!"); 332 return *levelIndexSets_[level]; 333 } 334 335 /** \brief Access to the LeafIndexSet */ leafIndexSet() const336 const typename Traits::LeafIndexSet& leafIndexSet() const 337 { 338 return leafIndexSet_; 339 } 340 341 /** @name Grid Refinement Methods */ 342 /*@{*/ 343 344 /** \brief Mark element for refinement 345 \param refCount <ul> 346 <li> 1: mark for red refinement </li> 347 <li> -1: mark for coarsening </li> 348 <li> 0: delete a possible refinement mark </li> 349 </ul> 350 \param e Element to be marked 351 \return <ul> 352 <li> true, if element was marked </li> 353 <li> false, if nothing changed </li> 354 </ul> 355 */ 356 bool mark(int refCount, const typename Traits::template Codim<0>::Entity & e ); 357 358 /** \brief Mark method accepting a UG refinement rule 359 360 \param e element to be marked for refinement 361 \param rule One of the UG refinement rules 362 \param side If rule==UG::%D2::%BLUE (one quadrilateral is split into two rectangles) 363 you can choose the orientation of the cut by setting side==0 or side==1 364 365 The available values for RefinementRule are: (see the RefinementRule enum in ug/gm/gm.h) 366 <h3>2D</h3> 367 368 - NO_REFINEMENT 369 - COPY 370 - RED 371 - BLUE 372 - COARSE 373 - BISECTION_1 374 - BISECTION_2_Q 375 - BISECTION_2_T1 376 - BISECTION_2_T2 377 - BISECTION_3 378 379 <h3>3D</h3> 380 381 - NO_REFINEMENT 382 - COPY 383 - RED 384 - COARSE 385 386 - TETRA_RED_HEX 387 388 - PRISM_BISECT_1_2 389 - PRISM_QUADSECT 390 - PRISM_BISECT_HEX0 391 - PRISM_BISECT_HEX1 392 - PRISM_BISECT_HEX2 393 - PRISM_ROTATE_LEFT 394 - PRISM_ROTATE_RGHT 395 - PRISM_QUADSECT_HEXPRI0 396 - PRISM_RED_HEX 397 - PRISM_BISECT_0_1 398 - PRISM_BISECT_0_2 399 - PRISM_BISECT_0_3 400 401 - HEX_BISECT_0_1 402 - HEX_BISECT_0_2 403 - HEX_BISECT_0_3 404 - HEX_TRISECT_0 405 - HEX_TRISECT_5 406 - HEX_QUADSECT_0 407 - HEX_QUADSECT_1 408 - HEX_QUADSECT_2 409 - HEX_BISECT_HEXPRI0 410 - HEX_BISECT_HEXPRI1 411 412 */ 413 bool mark(const typename Traits::template Codim<0>::Entity & e, 414 typename UG_NS<dim>::RefinementRule rule, 415 int side=0); 416 417 /** \brief Query whether element is marked for refinement */ 418 int getMark(const typename Traits::template Codim<0>::Entity& e) const; 419 420 /** \brief returns true, if some elements might be coarsend during grid 421 adaption, here always returns true */ 422 bool preAdapt(); 423 424 //! Triggers the grid refinement process 425 bool adapt(); 426 427 /** \brief Clean up refinement markers */ 428 void postAdapt(); 429 /*@}*/ 430 431 /** \brief Distributes the grid and some data over the available nodes in a distributed machine 432 433 \tparam DataHandle works like the data handle for the communicate 434 methods. 435 436 \return True, if grid has changed, false otherwise 437 */ 438 template<class DataHandle> loadBalance(DataHandle & dataHandle)439 bool loadBalance (DataHandle& dataHandle) 440 { 441 #ifdef ModelP 442 // gather element data 443 if (dataHandle.contains(dim, 0)) 444 UGLBGatherScatter::template gather<0>(this->leafGridView(), dataHandle); 445 446 // gather node data 447 if (dataHandle.contains(dim,dim)) 448 UGLBGatherScatter::template gather<dim>(this->leafGridView(), dataHandle); 449 #endif 450 451 // the load balancing step now also attaches 452 // the data to the entities and distributes it 453 loadBalance(); 454 455 #ifdef ModelP 456 // scatter element data 457 if (dataHandle.contains(dim, 0)) 458 UGLBGatherScatter::template scatter<0>(this->leafGridView(), dataHandle); 459 460 // scatter node data 461 if (dataHandle.contains(dim,dim)) 462 UGLBGatherScatter::template scatter<dim>(this->leafGridView(), dataHandle); 463 #endif 464 465 return true; 466 } 467 468 /** \brief Distributes this grid over the available nodes in a distributed machine 469 470 \bug The return value is always 'true' 471 472 \param minlevel The coarsest grid level that gets distributed 473 */ 474 bool loadBalance(int minlevel=0); 475 476 /** \brief Distribute this grid over a distributed machine 477 * 478 * \param[in] targetProcessors For each leaf element the rank of the process the element shall be sent to 479 * \param[in] fromLevel The lowest level that gets redistributed (set to 0 when in doubt) 480 * 481 * This method allows to (re-)distribute the grid controlled by an external grid repartitioning library. 482 * You need to get that library to assign a target rank to each interior element in the leaf grid. With this 483 * information in a std::vector, call this method, and UG will do the actual repartitioning. 484 * Each leaf element will be sent to the assigned target rank. For all other elements we look at 485 * where there children are being sent to. The parent is then sent to where most of its children are 486 * (Familienzusammenfuehrung). 487 * 488 * The size of the input array targetProcessors is expected to be equal to the number of elements in 489 * the 'all'-partition, i.e., the number Interior elements plus the number of Ghost elements. 490 * To get the array entry corresponding to an Interior element, a MultipleCodimMultipleGeomTypeMapper 491 * with layout class MCMGElementLayout is used. Array entries corresponding to Ghost elements are ignored. 492 * 493 * In some cases you may also want to leave the lowest levels on one process, to have them all together 494 * for multigrid coarse grid corrections. In that case, use the fromLevel parameter with a value other 495 * than zero, to redistribute only elements above a certain level. 496 * 497 * The fromLevel argument is also needed to allow the compiler to distinguish this method from 498 * the loadBalance method with a single template DataHandle argument. 499 * 500 * \note In theory you can assign a target rank to any element on any level, and UG will magically transfer 501 * the element to that rank and make everything come out right. This is not supported by the UGGrid interface, 502 * because I didn't see a use case for it. If you do need it please ask on the Dune mailing list. 503 * 504 * \return true 505 */ 506 bool loadBalance(const std::vector<Rank>& targetProcessors, unsigned int fromLevel); 507 508 /** \brief Distributes the grid over the processes of a parallel machine, and sends data along with it 509 * 510 * \param[in] targetProcessors For each leaf element the rank of the process the element shall be sent to 511 * \param[in] fromLevel The lowest level that gets redistributed (set to 0 when in doubt) 512 * \param[in,out] dataHandle A data handle object that does the gathering and scattering of data 513 * \tparam DataHandle works like the data handle for the communicate methods. 514 * 515 * \return true 516 */ 517 template<class DataHandle> loadBalance(const std::vector<Rank> & targetProcessors,unsigned int fromLevel,DataHandle & dataHandle)518 bool loadBalance (const std::vector<Rank>& targetProcessors, unsigned int fromLevel, DataHandle& dataHandle) 519 { 520 #ifdef ModelP 521 // gather element data 522 if (dataHandle.contains(dim, 0)) 523 UGLBGatherScatter::template gather<0>(this->leafGridView(), dataHandle); 524 525 // gather node data 526 if (dataHandle.contains(dim,dim)) 527 UGLBGatherScatter::template gather<dim>(this->leafGridView(), dataHandle); 528 #endif 529 530 // the load balancing step now also attaches 531 // the data to the entities and distributes it 532 loadBalance(targetProcessors,fromLevel); 533 534 #ifdef ModelP 535 // scatter element data 536 if (dataHandle.contains(dim, 0)) 537 UGLBGatherScatter::template scatter<0>(this->leafGridView(), dataHandle); 538 539 // scatter node data 540 if (dataHandle.contains(dim,dim)) 541 UGLBGatherScatter::template scatter<dim>(this->leafGridView(), dataHandle); 542 #endif 543 544 return true; 545 } 546 547 /** the collective communication */ comm() const548 const UGCollectiveCommunication& comm () const 549 { 550 return ccobj_; 551 } 552 553 protected: 554 #ifdef ModelP 555 template <int codim, class GridView, class DataHandle> communicateUG_(const GridView & gv,int level,DataHandle & dataHandle,InterfaceType iftype,CommunicationDirection dir) const556 void communicateUG_(const GridView& gv, int level, 557 DataHandle &dataHandle, 558 InterfaceType iftype, 559 CommunicationDirection dir) const 560 { 561 typename UG_NS<dim>::DDD_IF_DIR ugIfDir; 562 // Translate the communication direction from Dune-Speak to UG-Speak 563 if (dir==ForwardCommunication) 564 ugIfDir = UG_NS<dim>::IF_FORWARD(); 565 else 566 ugIfDir = UG_NS<dim>::IF_BACKWARD(); 567 568 typedef UGMessageBuffer<DataHandle,dim,codim> UGMsgBuf; 569 UGMsgBuf::duneDataHandle_ = &dataHandle; 570 571 UGMsgBuf::level = level; 572 573 std::vector<typename UG_NS<dim>::DDD_IF> ugIfs = findDDDInterfaces(iftype, codim); 574 575 unsigned bufSize = UGMsgBuf::ugBufferSize(gv); 576 if (!bufSize) 577 return; // we don't need to communicate if we don't have any data! 578 UGMsgBuf::grid_ = this; 579 for (unsigned i=0; i < ugIfs.size(); ++i) 580 UG_NS<dim>::DDD_IFOneway(multigrid_->dddContext(), 581 ugIfs[i], 582 ugIfDir, 583 bufSize, 584 &UGMsgBuf::ugGather_, 585 &UGMsgBuf::ugScatter_); 586 } 587 588 /** \brief Translate Dune communication interface to UG communication interface */ 589 std::vector<typename UG_NS<dim>::DDD_IF> findDDDInterfaces(InterfaceType iftype, 590 int codim) const; 591 #endif 592 public: 593 // ********************************************************** 594 // End of Interface Methods 595 // ********************************************************** 596 597 /** \brief Rudimentary substitute for a hierarchic iterator on faces 598 \param e, elementSide Grid face specified by an element and one of its sides 599 \param maxl The finest level that should be traversed by the iterator 600 \param[out] childElements For each subface: element index, elementSide, and level 601 \param[out] childElementSides Indices for transformation because Dune numbers the 602 faces of several elements differently than UG 603 */ 604 void getChildrenOfSubface(const typename Traits::template Codim<0>::Entity & e, 605 int elementSide, 606 int maxl, 607 std::vector<typename Traits::template Codim<0>::Entity>& childElements, 608 std::vector<unsigned char>& childElementSides) const; 609 610 /** \brief The different forms of grid refinement that UG supports */ 611 enum RefinementType { 612 /** \brief New level consists only of the refined elements and the closure*/ 613 LOCAL, 614 /** \brief New level consists of the refined elements and the unrefined ones, too */ 615 COPY 616 }; 617 618 /** \brief Decide whether to add a green closure to locally refined grid sections or not */ 619 enum ClosureType { 620 /** \brief Standard red/green refinement */ 621 GREEN, 622 /** \brief No closure, results in nonconforming meshes */ 623 NONE 624 }; 625 626 /** \brief Sets the type of grid refinement */ setRefinementType(RefinementType type)627 void setRefinementType(RefinementType type) { 628 refinementType_ = type; 629 } 630 631 /** \brief Sets the type of grid refinement closure */ setClosureType(ClosureType type)632 void setClosureType(ClosureType type) { 633 closureType_ = type; 634 } 635 636 /** \brief Sets a vertex to a new position 637 638 Changing a vertex' position changes its position on all grid levels!*/ 639 void setPosition(const typename Traits::template Codim<dim>::Entity& e, 640 const FieldVector<double, dim>& pos); 641 642 /** \brief Does uniform refinement 643 * 644 * \param n Number of uniform refinement steps 645 */ 646 void globalRefine(int n); 647 648 /** \brief Save entire grid hierarchy to disk 649 650 Test implementation -- not working! 651 */ 652 void saveState(const std::string& filename) const; 653 654 /** \brief Read entire grid hierarchy from disk 655 656 Test implementation -- not working! 657 */ 658 void loadState(const std::string& filename); 659 660 private: 661 /** \brief UG multigrid, which contains the actual grid hierarchy structure */ 662 typename UG_NS<dim>::MultiGrid* multigrid_; 663 664 /** \brief The collective communication object. */ 665 UGCollectiveCommunication ccobj_; 666 667 /** \brief Recomputes entity indices after the grid was changed 668 \param setLevelZero If this is false, level indices of the level 0 are not touched 669 \param nodePermutation Permutation array for the vertex level 0 indices. If this is NULL, 670 the identity is used. 671 */ 672 void setIndices(bool setLevelZero, 673 std::vector<unsigned int>* nodePermutation); 674 675 // Each UGGrid object has a unique name to identify it in the 676 // UG environment structure 677 std::string name_; 678 679 // Our set of level indices 680 std::vector<std::shared_ptr<UGGridLevelIndexSet<const UGGrid<dim> > > > levelIndexSets_; 681 682 UGGridLeafIndexSet<const UGGrid<dim> > leafIndexSet_; 683 684 // One id set implementation 685 // Used for both the local and the global UGGrid id sets 686 UGGridIdSet<const UGGrid<dim> > idSet_; 687 688 //! The type of grid refinement currently in use 689 RefinementType refinementType_; 690 691 //! The type of grid refinement closure currently in use 692 ClosureType closureType_; 693 694 /** \brief Number of UGGrids currently in use. 695 * 696 * This counts the number of UGGrids currently instantiated. All 697 * constructors of UGGrid look at this variable. If it is zero, they 698 * initialize UG before proceeding. Destructors use the same mechanism 699 * to safely shut down UG after deleting the last UGGrid object. 700 */ 701 static int numOfUGGrids; 702 703 /** \brief Remember whether some element has been marked for refinement 704 ever since the last call to adapt(). 705 706 This is here to implement the return value of adapt(). 707 */ 708 bool someElementHasBeenMarkedForRefinement_; 709 710 /** \brief Remember whether some element has been marked for coarsening 711 ever since the last call to adapt(). 712 713 This is here to implement the return value of preAdapt(). 714 */ 715 bool someElementHasBeenMarkedForCoarsening_; 716 717 /** \brief The classes implementing the geometry of the boundary segments, if requested */ 718 std::vector<std::shared_ptr<BoundarySegment<dim> > > boundarySegments_; 719 720 /** \brief Overall number of coarse grid boundary segments. 721 722 This includes the number of linear segments. 723 Hence numBoundarySegments_ >= boundarySegments_.size() (greater than or equal) 724 */ 725 unsigned int numBoundarySegments_; 726 727 }; // end Class UGGrid 728 729 namespace Capabilities 730 { 731 /** \struct hasEntity 732 \ingroup UGGrid 733 */ 734 735 /** \struct hasBackupRestoreFacilities 736 \ingroup UGGrid 737 */ 738 739 /** \struct IsUnstructured 740 \ingroup UGGrid 741 */ 742 743 /** \brief UGGrid has entities of all codimensions 744 \ingroup UGGrid 745 */ 746 template<int dim, int codim> 747 struct hasEntity< UGGrid<dim>, codim> 748 { 749 static const bool v = true; 750 }; 751 752 /** 753 * \brief Set default for hasEntityIterator to false 754 * UGGrid can currently only iterate over elements and vertices 755 * \ingroup UGGrid 756 **/ 757 template<int dim, int codim> 758 struct hasEntityIterator<UGGrid<dim>, codim> 759 { 760 static const bool v = false; 761 }; 762 763 /** 764 * \brief UGGrid can iterate over codim=0 entities (elements) 765 * \ingroup UGGrid 766 **/ 767 template<int dim> 768 struct hasEntityIterator<UGGrid<dim>, 0> 769 { 770 static const bool v = true; 771 }; 772 773 /** 774 * \brief UGGrid can iterate over codim=dim entities (vertices) 775 * \ingroup UGGrid 776 **/ 777 template<int dim> 778 struct hasEntityIterator<UGGrid<dim>, dim> 779 { 780 static const bool v = true; 781 }; 782 783 /** \brief UGGrid can communicate on entities of all (existing) codimensions 784 * \ingroup UGGrid 785 */ 786 template<int dim, int codim> 787 struct canCommunicate<UGGrid<dim>, codim> 788 { 789 static const bool v = (codim>=0 && codim<=dim); 790 }; 791 792 /** \brief UGGrid is levelwise conforming 793 \ingroup UGGrid 794 */ 795 template<int dim> 796 struct isLevelwiseConforming< UGGrid<dim> > 797 { 798 static const bool v = true; 799 }; 800 801 /** \brief UGGrid may not be leafwise conforming 802 \ingroup UGGrid 803 */ 804 template<int dim> 805 struct isLeafwiseConforming< UGGrid<dim> > 806 { 807 static const bool v = false; 808 }; 809 810 } 811 812 } // namespace Dune 813 814 #endif // HAVE_UG || DOXYGEN 815 #endif // DUNE_UGGRID_HH 816