dune-grid  2.7.0
yaspgrid.hh
Go to the documentation of this file.
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_GRID_YASPGRID_HH
4 #define DUNE_GRID_YASPGRID_HH
5 
6 #include <iostream>
7 #include <vector>
8 #include <algorithm>
9 #include <stack>
10 #include <type_traits>
11 
12 // either include stdint.h or provide fallback for uint8_t
13 #if HAVE_STDINT_H
14 #include <stdint.h>
15 #else
16 typedef unsigned char uint8_t;
17 #endif
18 
20 #include <dune/grid/common/grid.hh> // the grid base classes
21 #include <dune/grid/common/capabilities.hh> // the capabilities
22 #include <dune/common/hybridutilities.hh>
23 #include <dune/common/power.hh>
24 #include <dune/common/bigunsignedint.hh>
25 #include <dune/common/typetraits.hh>
26 #include <dune/common/reservedvector.hh>
27 #include <dune/common/parallel/communication.hh>
28 #include <dune/common/parallel/mpihelper.hh>
29 #include <dune/common/deprecated.hh>
30 #include <dune/geometry/axisalignedcubegeometry.hh>
31 #include <dune/geometry/type.hh>
34 
35 
36 #if HAVE_MPI
37 #include <dune/common/parallel/mpicommunication.hh>
38 #endif
39 
47 namespace Dune {
48 
49  /* some sizes for building global ids
50  */
51  const int yaspgrid_dim_bits = 24; // bits for encoding each dimension
52  const int yaspgrid_level_bits = 5; // bits for encoding level number
53 
54 
55  //************************************************************************
56  // forward declaration of templates
57 
58  template<int dim, class Coordinates> class YaspGrid;
59  template<int mydim, int cdim, class GridImp> class YaspGeometry;
60  template<int codim, int dim, class GridImp> class YaspEntity;
61  template<int codim, class GridImp> class YaspEntitySeed;
62  template<int codim, PartitionIteratorType pitype, class GridImp> class YaspLevelIterator;
63  template<class GridImp> class YaspIntersectionIterator;
64  template<class GridImp> class YaspIntersection;
65  template<class GridImp> class YaspHierarchicIterator;
66  template<class GridImp, bool isLeafIndexSet> class YaspIndexSet;
67  template<class GridImp> class YaspGlobalIdSet;
68  template<class GridImp> class YaspPersistentContainerIndex;
69 
70 } // namespace Dune
71 
85 
86 namespace Dune {
87 
88  template<int dim, class Coordinates>
90  {
91 #if HAVE_MPI
92  typedef CollectiveCommunication<MPI_Comm> CCType;
93 #else
94  typedef CollectiveCommunication<No_Comm> CCType;
95 #endif
96 
97  typedef GridTraits<dim, // dimension of the grid
98  dim, // dimension of the world space
101  YaspLevelIterator, // type used for the level iterator
102  YaspIntersection, // leaf intersection
103  YaspIntersection, // level intersection
104  YaspIntersectionIterator, // leaf intersection iter
105  YaspIntersectionIterator, // level intersection iter
107  YaspLevelIterator, // type used for the leaf(!) iterator
108  YaspIndexSet< const YaspGrid< dim, Coordinates >, false >, // level index set
109  YaspIndexSet< const YaspGrid< dim, Coordinates >, true >, // leaf index set
111  bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
113  bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
114  CCType,
118  };
119 
120 #ifndef DOXYGEN
121  template<int dim, int codim>
122  struct YaspCommunicateMeta {
123  template<class G, class DataHandle>
124  static void comm (const G& g, DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level)
125  {
126  if (data.contains(dim,codim))
127  {
128  g.template communicateCodim<DataHandle,codim>(data,iftype,dir,level);
129  }
130  YaspCommunicateMeta<dim,codim-1>::comm(g,data,iftype,dir,level);
131  }
132  };
133 
134  template<int dim>
135  struct YaspCommunicateMeta<dim,0> {
136  template<class G, class DataHandle>
137  static void comm (const G& g, DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level)
138  {
139  if (data.contains(dim,0))
140  g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
141  }
142  };
143 #endif
144 
145  //************************************************************************
162  template<int dim, class Coordinates = EquidistantCoordinates<double, dim> >
163  class YaspGrid
164  : public GridDefaultImplementation<dim,dim,typename Coordinates::ctype,YaspGridFamily<dim, Coordinates> >
165  {
166 
167  template<int, PartitionIteratorType, typename>
168  friend class YaspLevelIterator;
169 
170  template<typename>
172 
173  protected:
174 
176 
177  public:
179  typedef typename Coordinates::ctype ctype;
180 #if HAVE_MPI
181  typedef CollectiveCommunication<MPI_Comm> CollectiveCommunicationType;
182 #else
183  typedef CollectiveCommunication<No_Comm> CollectiveCommunicationType;
184 #endif
185 
186 #ifndef DOXYGEN
187  typedef typename Dune::YGrid<Coordinates> YGrid;
189 
192  struct YGridLevel {
193 
195  int level() const
196  {
197  return level_;
198  }
199 
200  Coordinates coords;
201 
202  std::array<YGrid, dim+1> overlapfront;
203  std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlapfront_data;
204  std::array<YGrid, dim+1> overlap;
205  std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlap_data;
206  std::array<YGrid, dim+1> interiorborder;
207  std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interiorborder_data;
208  std::array<YGrid, dim+1> interior;
209  std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interior_data;
210 
211  std::array<YGridList<Coordinates>,dim+1> send_overlapfront_overlapfront;
212  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlapfront_overlapfront_data;
213  std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlapfront;
214  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlapfront_data;
215 
216  std::array<YGridList<Coordinates>,dim+1> send_overlap_overlapfront;
217  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlap_overlapfront_data;
218  std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlap;
219  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlap_data;
220 
221  std::array<YGridList<Coordinates>,dim+1> send_interiorborder_interiorborder;
222  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_interiorborder_data;
223  std::array<YGridList<Coordinates>,dim+1> recv_interiorborder_interiorborder;
224  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_interiorborder_interiorborder_data;
225 
226  std::array<YGridList<Coordinates>,dim+1> send_interiorborder_overlapfront;
227  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_overlapfront_data;
228  std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_interiorborder;
229  std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_interiorborder_data;
230 
231  // general
232  YaspGrid<dim,Coordinates>* mg; // each grid level knows its multigrid
233  int overlapSize; // in mesh cells on this level
234  bool keepOverlap;
235 
237  int level_;
238  };
239 
241  typedef std::array<int, dim> iTupel;
242  typedef FieldVector<ctype, dim> fTupel;
243 
244  // communication tag used by multigrid
245  enum { tag = 17 };
246 #endif
247 
250  {
251  return _torus;
252  }
253 
255  int globalSize(int i) const
256  {
257  return levelSize(maxLevel(),i);
258  }
259 
261  iTupel globalSize() const
262  {
263  return levelSize(maxLevel());
264  }
265 
267  int levelSize(int l, int i) const
268  {
269  return _coarseSize[i] * (1 << l);
270  }
271 
273  iTupel levelSize(int l) const
274  {
275  iTupel s;
276  for (int i=0; i<dim; ++i)
277  s[i] = levelSize(l,i);
278  return s;
279  }
280 
282  bool isPeriodic(int i) const
283  {
284  return _periodic[i];
285  }
286 
287  bool getRefineOption() const
288  {
289  return keep_ovlp;
290  }
291 
293  typedef typename ReservedVector<YGridLevel,32>::const_iterator YGridLevelIterator;
294 
297  {
298  return YGridLevelIterator(_levels,0);
299  }
300 
302  YGridLevelIterator begin (int i) const
303  {
304  if (i<0 || i>maxLevel())
305  DUNE_THROW(GridError, "level not existing");
306  return YGridLevelIterator(_levels,i);
307  }
308 
311  {
312  return YGridLevelIterator(_levels,maxLevel()+1);
313  }
314 
315  // static method to create the default load balance strategy
317  {
318  static YLoadBalanceDefault<dim> lb;
319  return & lb;
320  }
321 
322  protected:
330  void makelevel (const Coordinates& coords, std::bitset<dim> periodic, iTupel o_interior, int overlap)
331  {
332  YGridLevel& g = _levels.back();
333  g.overlapSize = overlap;
334  g.mg = this;
335  g.level_ = maxLevel();
336  g.coords = coords;
337  g.keepOverlap = keep_ovlp;
338 
339  // set the inserting positions in the corresponding arrays of YGridLevelStructure
340  typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlapfront_it = g.overlapfront_data.begin();
341  typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlap_it = g.overlap_data.begin();
342  typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interiorborder_it = g.interiorborder_data.begin();
343  typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interior_it = g.interior_data.begin();
344 
345  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
346  send_overlapfront_overlapfront_it = g.send_overlapfront_overlapfront_data.begin();
347  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
348  recv_overlapfront_overlapfront_it = g.recv_overlapfront_overlapfront_data.begin();
349 
350  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
351  send_overlap_overlapfront_it = g.send_overlap_overlapfront_data.begin();
352  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
353  recv_overlapfront_overlap_it = g.recv_overlapfront_overlap_data.begin();
354 
355  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
356  send_interiorborder_interiorborder_it = g.send_interiorborder_interiorborder_data.begin();
357  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
358  recv_interiorborder_interiorborder_it = g.recv_interiorborder_interiorborder_data.begin();
359 
360  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
361  send_interiorborder_overlapfront_it = g.send_interiorborder_overlapfront_data.begin();
362  typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
363  recv_overlapfront_interiorborder_it = g.recv_overlapfront_interiorborder_data.begin();
364 
365  // have a null array for constructor calls around
366  std::array<int,dim> n;
367  std::fill(n.begin(), n.end(), 0);
368 
369  // determine origin of the grid with overlap and store whether an overlap area exists in direction i.
370  std::bitset<dim> ovlp_low(0ULL);
371  std::bitset<dim> ovlp_up(0ULL);
372 
373  iTupel o_overlap;
374  iTupel s_overlap;
375 
376  // determine at where we have overlap and how big the size of the overlap partition is
377  for (int i=0; i<dim; i++)
378  {
379  // the coordinate container has been contructed to hold the entire grid on
380  // this processor, including overlap. this is the element size.
381  s_overlap[i] = g.coords.size(i);
382 
383  //in the periodic case there is always overlap
384  if (periodic[i])
385  {
386  o_overlap[i] = o_interior[i]-overlap;
387  ovlp_low[i] = true;
388  ovlp_up[i] = true;
389  }
390  else
391  {
392  //check lower boundary
393  if (o_interior[i] - overlap < 0)
394  o_overlap[i] = 0;
395  else
396  {
397  o_overlap[i] = o_interior[i] - overlap;
398  ovlp_low[i] = true;
399  }
400 
401  //check upper boundary
402  if (o_overlap[i] + g.coords.size(i) < globalSize(i))
403  ovlp_up[i] = true;
404  }
405  }
406 
407  for (unsigned int codim = 0; codim < dim + 1; codim++)
408  {
409  // set the begin iterator for the corresponding ygrids
410  g.overlapfront[codim].setBegin(overlapfront_it);
411  g.overlap[codim].setBegin(overlap_it);
412  g.interiorborder[codim].setBegin(interiorborder_it);
413  g.interior[codim].setBegin(interior_it);
414  g.send_overlapfront_overlapfront[codim].setBegin(send_overlapfront_overlapfront_it);
415  g.recv_overlapfront_overlapfront[codim].setBegin(recv_overlapfront_overlapfront_it);
416  g.send_overlap_overlapfront[codim].setBegin(send_overlap_overlapfront_it);
417  g.recv_overlapfront_overlap[codim].setBegin(recv_overlapfront_overlap_it);
418  g.send_interiorborder_interiorborder[codim].setBegin(send_interiorborder_interiorborder_it);
419  g.recv_interiorborder_interiorborder[codim].setBegin(recv_interiorborder_interiorborder_it);
420  g.send_interiorborder_overlapfront[codim].setBegin(send_interiorborder_overlapfront_it);
421  g.recv_overlapfront_interiorborder[codim].setBegin(recv_overlapfront_interiorborder_it);
422 
423  // find all combinations of unit vectors that span entities of the given codimension
424  for (unsigned int index = 0; index < (1<<dim); index++)
425  {
426  // check whether the given shift is of our codimension
427  std::bitset<dim> r(index);
428  if (r.count() != dim-codim)
429  continue;
430 
431  // get an origin and a size array for subsequent modification
432  std::array<int,dim> origin(o_overlap);
433  std::array<int,dim> size(s_overlap);
434 
435  // build overlapfront
436  // we have to extend the element size by one in all directions without shift.
437  for (int i=0; i<dim; i++)
438  if (!r[i])
439  size[i]++;
440  *overlapfront_it = YGridComponent<Coordinates>(origin, r, &g.coords, size, n, size);
441 
442  // build overlap
443  for (int i=0; i<dim; i++)
444  {
445  if (!r[i])
446  {
447  if (ovlp_low[i])
448  {
449  origin[i]++;
450  size[i]--;
451  }
452  if (ovlp_up[i])
453  size[i]--;
454  }
455  }
456  *overlap_it = YGridComponent<Coordinates>(origin,size,*overlapfront_it);
457 
458  // build interiorborder
459  for (int i=0; i<dim; i++)
460  {
461  if (ovlp_low[i])
462  {
463  origin[i] += overlap;
464  size[i] -= overlap;
465  if (!r[i])
466  {
467  origin[i]--;
468  size[i]++;
469  }
470  }
471  if (ovlp_up[i])
472  {
473  size[i] -= overlap;
474  if (!r[i])
475  size[i]++;
476  }
477  }
478  *interiorborder_it = YGridComponent<Coordinates>(origin,size,*overlapfront_it);
479 
480  // build interior
481  for (int i=0; i<dim; i++)
482  {
483  if (!r[i])
484  {
485  if (ovlp_low[i])
486  {
487  origin[i]++;
488  size[i]--;
489  }
490  if (ovlp_up[i])
491  size[i]--;
492  }
493  }
494  *interior_it = YGridComponent<Coordinates>(origin, size, *overlapfront_it);
495 
496  intersections(*overlapfront_it,*overlapfront_it,*send_overlapfront_overlapfront_it, *recv_overlapfront_overlapfront_it);
497  intersections(*overlap_it,*overlapfront_it,*send_overlap_overlapfront_it, *recv_overlapfront_overlap_it);
498  intersections(*interiorborder_it,*interiorborder_it,*send_interiorborder_interiorborder_it,*recv_interiorborder_interiorborder_it);
499  intersections(*interiorborder_it,*overlapfront_it,*send_interiorborder_overlapfront_it,*recv_overlapfront_interiorborder_it);
500 
501  // advance all iterators pointing to the next insertion point
502  ++overlapfront_it;
503  ++overlap_it;
504  ++interiorborder_it;
505  ++interior_it;
506  ++send_overlapfront_overlapfront_it;
507  ++recv_overlapfront_overlapfront_it;
508  ++send_overlap_overlapfront_it;
509  ++recv_overlapfront_overlap_it;
510  ++send_interiorborder_interiorborder_it;
511  ++recv_interiorborder_interiorborder_it;
512  ++send_interiorborder_overlapfront_it;
513  ++recv_overlapfront_interiorborder_it;
514  }
515 
516  // set end iterators in the corresonding ygrids
517  g.overlapfront[codim].finalize(overlapfront_it);
518  g.overlap[codim].finalize(overlap_it);
519  g.interiorborder[codim].finalize(interiorborder_it);
520  g.interior[codim].finalize(interior_it);
521  g.send_overlapfront_overlapfront[codim].finalize(send_overlapfront_overlapfront_it,g.overlapfront[codim]);
522  g.recv_overlapfront_overlapfront[codim].finalize(recv_overlapfront_overlapfront_it,g.overlapfront[codim]);
523  g.send_overlap_overlapfront[codim].finalize(send_overlap_overlapfront_it,g.overlapfront[codim]);
524  g.recv_overlapfront_overlap[codim].finalize(recv_overlapfront_overlap_it,g.overlapfront[codim]);
525  g.send_interiorborder_interiorborder[codim].finalize(send_interiorborder_interiorborder_it,g.overlapfront[codim]);
526  g.recv_interiorborder_interiorborder[codim].finalize(recv_interiorborder_interiorborder_it,g.overlapfront[codim]);
527  g.send_interiorborder_overlapfront[codim].finalize(send_interiorborder_overlapfront_it,g.overlapfront[codim]);
528  g.recv_overlapfront_interiorborder[codim].finalize(recv_overlapfront_interiorborder_it,g.overlapfront[codim]);
529  }
530  }
531 
532 #ifndef DOXYGEN
533 
541  struct mpifriendly_ygrid {
542  mpifriendly_ygrid ()
543  {
544  std::fill(origin.begin(), origin.end(), 0);
545  std::fill(size.begin(), size.end(), 0);
546  }
547  mpifriendly_ygrid (const YGridComponent<Coordinates>& grid)
548  : origin(grid.origin()), size(grid.size())
549  {}
550  iTupel origin;
551  iTupel size;
552  };
553 #endif
554 
564  std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
565  {
566  iTupel size = globalSize();
567 
568  // the exchange buffers
569  std::vector<YGridComponent<Coordinates> > send_recvgrid(_torus.neighbors());
570  std::vector<YGridComponent<Coordinates> > recv_recvgrid(_torus.neighbors());
571  std::vector<YGridComponent<Coordinates> > send_sendgrid(_torus.neighbors());
572  std::vector<YGridComponent<Coordinates> > recv_sendgrid(_torus.neighbors());
573 
574  // new exchange buffers to send simple struct without virtual functions
575  std::vector<mpifriendly_ygrid> mpifriendly_send_recvgrid(_torus.neighbors());
576  std::vector<mpifriendly_ygrid> mpifriendly_recv_recvgrid(_torus.neighbors());
577  std::vector<mpifriendly_ygrid> mpifriendly_send_sendgrid(_torus.neighbors());
578  std::vector<mpifriendly_ygrid> mpifriendly_recv_sendgrid(_torus.neighbors());
579 
580  // fill send buffers; iterate over all neighboring processes
581  // non-periodic case is handled automatically because intersection will be zero
582  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
583  {
584  // determine if we communicate with this neighbor (and what)
585  bool skip = false;
586  iTupel coord = _torus.coord(); // my coordinates
587  iTupel delta = i.delta(); // delta to neighbor
588  iTupel nb = coord; // the neighbor
589  for (int k=0; k<dim; k++) nb[k] += delta[k];
590  iTupel v; // grid movement
591  std::fill(v.begin(), v.end(), 0);
592 
593  for (int k=0; k<dim; k++)
594  {
595  if (nb[k]<0)
596  {
597  if (_periodic[k])
598  v[k] += size[k];
599  else
600  skip = true;
601  }
602  if (nb[k]>=_torus.dims(k))
603  {
604  if (_periodic[k])
605  v[k] -= size[k];
606  else
607  skip = true;
608  }
609  // neither might be true, then v=0
610  }
611 
612  // store moved grids in send buffers
613  if (!skip)
614  {
615  send_sendgrid[i.index()] = sendgrid.move(v);
616  send_recvgrid[i.index()] = recvgrid.move(v);
617  }
618  else
619  {
620  send_sendgrid[i.index()] = YGridComponent<Coordinates>();
621  send_recvgrid[i.index()] = YGridComponent<Coordinates>();
622  }
623  }
624 
625  // issue send requests for sendgrid being sent to all neighbors
626  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
627  {
628  mpifriendly_send_sendgrid[i.index()] = mpifriendly_ygrid(send_sendgrid[i.index()]);
629  _torus.send(i.rank(), &mpifriendly_send_sendgrid[i.index()], sizeof(mpifriendly_ygrid));
630  }
631 
632  // issue recv requests for sendgrids of neighbors
633  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
634  _torus.recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()], sizeof(mpifriendly_ygrid));
635 
636  // exchange the sendgrids
637  _torus.exchange();
638 
639  // issue send requests for recvgrid being sent to all neighbors
640  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
641  {
642  mpifriendly_send_recvgrid[i.index()] = mpifriendly_ygrid(send_recvgrid[i.index()]);
643  _torus.send(i.rank(), &mpifriendly_send_recvgrid[i.index()], sizeof(mpifriendly_ygrid));
644  }
645 
646  // issue recv requests for recvgrid of neighbors
647  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
648  _torus.recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()], sizeof(mpifriendly_ygrid));
649 
650  // exchange the recvgrid
651  _torus.exchange();
652 
653  // process receive buffers and compute intersections
654  for (typename Torus<CollectiveCommunicationType,dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
655  {
656  // what must be sent to this neighbor
657  Intersection send_intersection;
658  mpifriendly_ygrid yg = mpifriendly_recv_recvgrid[i.index()];
659  recv_recvgrid[i.index()] = YGridComponent<Coordinates>(yg.origin,yg.size);
660  send_intersection.grid = sendgrid.intersection(recv_recvgrid[i.index()]);
661  send_intersection.rank = i.rank();
662  send_intersection.distance = i.distance();
663  if (!send_intersection.grid.empty()) sendlist.push_front(send_intersection);
664 
665  Intersection recv_intersection;
666  yg = mpifriendly_recv_sendgrid[i.index()];
667  recv_sendgrid[i.index()] = YGridComponent<Coordinates>(yg.origin,yg.size);
668  recv_intersection.grid = recvgrid.intersection(recv_sendgrid[i.index()]);
669  recv_intersection.rank = i.rank();
670  recv_intersection.distance = i.distance();
671  if(!recv_intersection.grid.empty()) recvlist.push_back(recv_intersection);
672  }
673  }
674 
675  protected:
676 
678 
679  void init()
680  {
681  indexsets.push_back( std::make_shared< YaspIndexSet<const YaspGrid<dim, Coordinates>, false > >(*this,0) );
683  }
684 
686  {
687  // sizes of local macro grid
688  std::array<int, dim> sides;
689  {
690  for (int i=0; i<dim; i++)
691  {
692  sides[i] =
693  ((begin()->overlap[0].dataBegin()->origin(i) == 0)+
694  (begin()->overlap[0].dataBegin()->origin(i) + begin()->overlap[0].dataBegin()->size(i)
695  == levelSize(0,i)));
696  }
697  }
698  nBSegments = 0;
699  for (int k=0; k<dim; k++)
700  {
701  int offset = 1;
702  for (int l=0; l<dim; l++)
703  {
704  if (l==k) continue;
705  offset *= begin()->overlap[0].dataBegin()->size(l);
706  }
707  nBSegments += sides[k]*offset;
708  }
709  }
710 
711  public:
712 
713  // define the persistent index type
714  typedef bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim> PersistentIndexType;
715 
718  // the Traits
720 
721  // need for friend declarations in entity
725 
733  YaspGrid (const Coordinates& coordinates,
734  std::bitset<dim> periodic = std::bitset<dim>(0ULL),
735  int overlap = 1,
738  : ccobj(comm)
739  , leafIndexSet_(*this)
740  , _periodic(periodic)
741  , _overlap(overlap)
742  , keep_ovlp(true)
743  , adaptRefCount(0)
744  , adaptActive(false)
745  {
746  _levels.resize(1);
747 
748  // Number of elements per coordinate direction on the coarsest level
749  for (std::size_t i=0; i<dim; i++)
750  _coarseSize[i] = coordinates.size(i);
751 
752  // Construct the communication torus
753  _torus = decltype(_torus)(comm,tag,_coarseSize,lb);
754 
755  iTupel o;
756  std::fill(o.begin(), o.end(), 0);
757  iTupel o_interior(o);
758  iTupel s_interior(_coarseSize);
759 
760  _torus.partition(_torus.rank(),o,_coarseSize,o_interior,s_interior);
761 
762  // Set domain size
763  if (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value
764  || std::is_same<Coordinates,EquidistantOffsetCoordinates<ctype,dim> >::value)
765  {
766  for (std::size_t i=0; i<dim; i++)
767  _L[i] = coordinates.size(i) * coordinates.meshsize(i,0);
768  }
769  if (std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >::value)
770  {
771  //determine sizes of vector to correctly construct torus structure and store for later size requests
772  for (int i=0; i<dim; i++)
773  _L[i] = coordinates.coordinate(i,_coarseSize[i]) - coordinates.coordinate(i,0);
774  }
775 
776 #if HAVE_MPI
777  // TODO: Settle on a single value for all coordinate types
778  int mysteryFactor = (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value) ? 1 : 2;
779 
780  // check whether the grid is large enough to be overlapping
781  for (int i=0; i<dim; i++)
782  {
783  // find out whether the grid is too small to
784  int toosmall = (s_interior[i] <= mysteryFactor * overlap) && // interior is very small
785  (periodic[i] || (s_interior[i] != _coarseSize[i])); // there is an overlap in that direction
786  // communicate the result to all those processes to have all processors error out if one process failed.
787  int global = 0;
788  MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR, comm);
789  if (global)
790  DUNE_THROW(Dune::GridError,"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
791  " Note that this also holds for DOFs on subdomain boundaries."
792  " Increase grid elements or decrease overlap accordingly.");
793  }
794 #endif // #if HAVE_MPI
795 
796  if (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value
797  || std::is_same<Coordinates,EquidistantOffsetCoordinates<ctype,dim> >::value)
798  {
799  iTupel s_overlap(s_interior);
800  for (int i=0; i<dim; i++)
801  {
802  if ((o_interior[i] - overlap > 0) || (periodic[i]))
803  s_overlap[i] += overlap;
804  if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
805  s_overlap[i] += overlap;
806  }
807 
808  FieldVector<ctype,dim> upperRightWithOverlap;
809  for (int i=0; i<dim; i++)
810  upperRightWithOverlap[i] = coordinates.coordinate(i,0) + coordinates.meshsize(i,0) * s_overlap[i];
811 
812  Hybrid::ifElse(std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >{}, [&](auto id)
813  {
814  // New coordinate object that additionally contains the overlap elements
815  EquidistantCoordinates<ctype,dim> coordinatesWithOverlap(upperRightWithOverlap,s_overlap);
816 
817  // add level (the this-> is needed to make g++-6 happy)
818  this->makelevel(id(coordinatesWithOverlap),periodic,o_interior,overlap);
819  });
820 
821  Hybrid::ifElse(std::is_same<Coordinates,EquidistantOffsetCoordinates<ctype,dim> >{}, [&](auto id)
822  {
823  Dune::FieldVector<ctype,dim> lowerleft;
824  for (int i=0; i<dim; i++)
825  lowerleft[i] = id(coordinates).origin(i);
826 
827  // New coordinate object that additionally contains the overlap elements
828  EquidistantOffsetCoordinates<ctype,dim> coordinatesWithOverlap(lowerleft,upperRightWithOverlap,s_overlap);
829 
830  // add level (the this-> is needed to make g++-6 happy)
831  this->makelevel(id(coordinatesWithOverlap),periodic,o_interior,overlap);
832  });
833  }
834 
835  Hybrid::ifElse(std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >{}, [&](auto id)
836  {
837  std::array<std::vector<ctype>,dim> newCoords;
838  std::array<int, dim> offset(o_interior);
839 
840  // find the relevant part of the coords vector for this processor and copy it to newCoords
841  for (int i=0; i<dim; ++i)
842  {
843  //define the coordinate range to be used
844  std::size_t begin = o_interior[i];
845  std::size_t end = begin + s_interior[i] + 1;
846 
847  // check whether we are not at the physical boundary. In that case overlap is a simple
848  // extension of the coordinate range to be used
849  if (o_interior[i] - overlap > 0)
850  {
851  begin = begin - overlap;
852  offset[i] -= overlap;
853  }
854  if (o_interior[i] + s_interior[i] + overlap < _coarseSize[i])
855  end = end + overlap;
856 
857  //copy the selected part in the new coord vector
858  newCoords[i].resize(end-begin);
859  auto newCoordsIt = newCoords[i].begin();
860  for (std::size_t j=begin; j<end; j++)
861  {
862  *newCoordsIt = coordinates.coordinate(i, j);
863  newCoordsIt++;
864  }
865 
866  // Check whether we are at the physical boundary and have a periodic grid.
867  // In this case the coordinate vector has to be tweaked manually.
868  if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
869  {
870  // we need to add the first <overlap> cells to the end of newcoords
871  for (int j=0; j<overlap; ++j)
872  newCoords[i].push_back(newCoords[i].back() - coordinates.coordinate(i,j) + coordinates.coordinate(i,j+1));
873  }
874 
875  if ((periodic[i]) && (o_interior[i] - overlap <= 0))
876  {
877  offset[i] -= overlap;
878 
879  // we need to add the last <overlap> cells to the begin of newcoords
880  std::size_t reverseCounter = coordinates.size(i);
881  for (int j=0; j<overlap; ++j, --reverseCounter)
882  newCoords[i].insert(newCoords[i].begin(), newCoords[i].front()
883  - coordinates.coordinate(i,reverseCounter) + coordinates.coordinate(i,reverseCounter-1));
884  }
885  }
886 
887  TensorProductCoordinates<ctype,dim> coordinatesWithOverlap(newCoords, offset);
888 
889  // add level (the this-> is needed to make g++-6 happy)
890  this->makelevel(id(coordinatesWithOverlap),periodic,o_interior,overlap);
891  });
892 
893  init();
894  }
895 
904  YaspGrid (Dune::FieldVector<ctype, dim> L,
905  std::array<int, dim> s,
906  std::bitset<dim> periodic = std::bitset<dim>(0ULL),
907  int overlap = 1,
910  : ccobj(comm), _torus(comm,tag,s,lb), leafIndexSet_(*this),
911  _L(L), _periodic(periodic), _coarseSize(s), _overlap(overlap),
912  keep_ovlp(true), adaptRefCount(0), adaptActive(false)
913  {
914  // check whether YaspGrid has been given the correct template parameter
915  static_assert(std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value,
916  "YaspGrid coordinate container template parameter and given constructor values do not match!");
917 
918  _levels.resize(1);
919 
920  iTupel o;
921  std::fill(o.begin(), o.end(), 0);
922  iTupel o_interior(o);
923  iTupel s_interior(s);
924 
925  _torus.partition(_torus.rank(),o,s,o_interior,s_interior);
926 
927 #if HAVE_MPI
928  // check whether the grid is large enough to be overlapping
929  for (int i=0; i<dim; i++)
930  {
931  // find out whether the grid is too small to
932  int toosmall = (s_interior[i] / 2 <= overlap) && // interior is very small
933  (periodic[i] || (s_interior[i] != s[i])); // there is an overlap in that direction
934  // communicate the result to all those processes to have all processors error out if one process failed.
935  int global = 0;
936  MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR, comm);
937  if (global)
938  DUNE_THROW(Dune::GridError,"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
939  " Note that this also holds for DOFs on subdomain boundaries."
940  " Increase grid elements or decrease overlap accordingly.");
941  }
942 #endif // #if HAVE_MPI
943 
944  iTupel s_overlap(s_interior);
945  for (int i=0; i<dim; i++)
946  {
947  if ((o_interior[i] - overlap > 0) || (periodic[i]))
948  s_overlap[i] += overlap;
949  if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
950  s_overlap[i] += overlap;
951  }
952 
953  FieldVector<ctype,dim> upperRightWithOverlap;
954 
955  for (int i=0; i<dim; i++)
956  upperRightWithOverlap[i] = (L[i] / s[i]) * s_overlap[i];
957 
958  // New coordinate object that additionally contains the overlap elements
959  EquidistantCoordinates<ctype,dim> cc(upperRightWithOverlap,s_overlap);
960 
961  // add level
962  makelevel(cc,periodic,o_interior,overlap);
963 
964  init();
965  }
966 
976  YaspGrid (Dune::FieldVector<ctype, dim> lowerleft,
977  Dune::FieldVector<ctype, dim> upperright,
978  std::array<int, dim> s,
979  std::bitset<dim> periodic = std::bitset<dim>(0ULL),
980  int overlap = 1,
983  : ccobj(comm), _torus(comm,tag,s,lb), leafIndexSet_(*this),
984  _L(upperright - lowerleft),
985  _periodic(periodic), _coarseSize(s), _overlap(overlap),
986  keep_ovlp(true), adaptRefCount(0), adaptActive(false)
987  {
988  // check whether YaspGrid has been given the correct template parameter
989  static_assert(std::is_same<Coordinates,EquidistantOffsetCoordinates<ctype,dim> >::value,
990  "YaspGrid coordinate container template parameter and given constructor values do not match!");
991 
992  _levels.resize(1);
993 
994  iTupel o;
995  std::fill(o.begin(), o.end(), 0);
996  iTupel o_interior(o);
997  iTupel s_interior(s);
998 
999  _torus.partition(_torus.rank(),o,s,o_interior,s_interior);
1000 
1001 #if HAVE_MPI
1002  // check whether the grid is large enough to be overlapping
1003  for (int i=0; i<dim; i++)
1004  {
1005  // find out whether the grid is too small to
1006  int toosmall = (s_interior[i] / 2 <= overlap) && // interior is very small
1007  (periodic[i] || (s_interior[i] != s[i])); // there is an overlap in that direction
1008  // communicate the result to all those processes to have all processors error out if one process failed.
1009  int global = 0;
1010  MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR, comm);
1011  if (global)
1012  DUNE_THROW(Dune::GridError,"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1013  " Note that this also holds for DOFs on subdomain boundaries."
1014  " Increase grid elements or decrease overlap accordingly.");
1015  }
1016 #endif // #if HAVE_MPI
1017 
1018  iTupel s_overlap(s_interior);
1019  for (int i=0; i<dim; i++)
1020  {
1021  if ((o_interior[i] - overlap > 0) || (periodic[i]))
1022  s_overlap[i] += overlap;
1023  if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
1024  s_overlap[i] += overlap;
1025  }
1026 
1027  FieldVector<ctype,dim> upperRightWithOverlap;
1028  for (int i=0; i<dim; i++)
1029  upperRightWithOverlap[i] = lowerleft[i]
1030  + s_overlap[i] * (upperright[i]-lowerleft[i]) / s[i];
1031 
1032  EquidistantOffsetCoordinates<ctype,dim> cc(lowerleft,upperRightWithOverlap,s_overlap);
1033 
1034  // add level
1035  makelevel(cc,periodic,o_interior,overlap);
1036 
1037  init();
1038  }
1039 
1047  YaspGrid (std::array<std::vector<ctype>, dim> coords,
1048  std::bitset<dim> periodic = std::bitset<dim>(0ULL),
1049  int overlap = 1,
1051  const YLoadBalance<dim>* lb = defaultLoadbalancer())
1052  : ccobj(comm), _torus(comm,tag,Dune::Yasp::sizeArray<dim>(coords),lb),
1053  leafIndexSet_(*this), _periodic(periodic), _overlap(overlap),
1054  keep_ovlp(true), adaptRefCount(0), adaptActive(false)
1055  {
1056  if (!Dune::Yasp::checkIfMonotonous(coords))
1057  DUNE_THROW(Dune::GridError,"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1058 
1059  // check whether YaspGrid has been given the correct template parameter
1060  static_assert(std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >::value,
1061  "YaspGrid coordinate container template parameter and given constructor values do not match!");
1062 
1063  _levels.resize(1);
1064 
1065  //determine sizes of vector to correctly construct torus structure and store for later size requests
1066  for (int i=0; i<dim; i++) {
1067  _coarseSize[i] = coords[i].size() - 1;
1068  _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
1069  }
1070 
1071  iTupel o;
1072  std::fill(o.begin(), o.end(), 0);
1073  iTupel o_interior(o);
1074  iTupel s_interior(_coarseSize);
1075 
1076  _torus.partition(_torus.rank(),o,_coarseSize,o_interior,s_interior);
1077 
1078 #if HAVE_MPI
1079  // check whether the grid is large enough to be overlapping
1080  for (int i=0; i<dim; i++)
1081  {
1082  // find out whether the grid is too small to
1083  int toosmall = (s_interior[i] / 2 <= overlap) && // interior is very small
1084  (periodic[i] || (s_interior[i] != _coarseSize[i])); // there is an overlap in that direction
1085  // communicate the result to all those processes to have all processors error out if one process failed.
1086  int global = 0;
1087  MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR, comm);
1088  if (global)
1089  DUNE_THROW(Dune::GridError,"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1090  " Note that this also holds for DOFs on subdomain boundaries."
1091  " Increase grid elements or decrease overlap accordingly.");
1092  }
1093 #endif // #if HAVE_MPI
1094 
1095 
1096  std::array<std::vector<ctype>,dim> newcoords;
1097  std::array<int, dim> offset(o_interior);
1098 
1099  // find the relevant part of the coords vector for this processor and copy it to newcoords
1100  for (int i=0; i<dim; ++i)
1101  {
1102  //define iterators on coords that specify the coordinate range to be used
1103  typename std::vector<ctype>::iterator begin = coords[i].begin() + o_interior[i];
1104  typename std::vector<ctype>::iterator end = begin + s_interior[i] + 1;
1105 
1106  // check whether we are not at the physical boundary. In that case overlap is a simple
1107  // extension of the coordinate range to be used
1108  if (o_interior[i] - overlap > 0)
1109  {
1110  begin = begin - overlap;
1111  offset[i] -= overlap;
1112  }
1113  if (o_interior[i] + s_interior[i] + overlap < _coarseSize[i])
1114  end = end + overlap;
1115 
1116  //copy the selected part in the new coord vector
1117  newcoords[i].resize(end-begin);
1118  std::copy(begin, end, newcoords[i].begin());
1119 
1120  // check whether we are at the physical boundary and a have a periodic grid.
1121  // In this case the coordinate vector has to be tweaked manually.
1122  if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
1123  {
1124  // we need to add the first <overlap> cells to the end of newcoords
1125  typename std::vector<ctype>::iterator it = coords[i].begin();
1126  for (int j=0; j<overlap; ++j)
1127  newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
1128  }
1129 
1130  if ((periodic[i]) && (o_interior[i] - overlap <= 0))
1131  {
1132  offset[i] -= overlap;
1133 
1134  // we need to add the last <overlap> cells to the begin of newcoords
1135  typename std::vector<ctype>::iterator it = coords[i].end() - 1;
1136  for (int j=0; j<overlap; ++j)
1137  newcoords[i].insert(newcoords[i].begin(), newcoords[i].front() - *it + *(--it));
1138  }
1139  }
1140 
1141  TensorProductCoordinates<ctype,dim> cc(newcoords, offset);
1142 
1143  // add level
1144  makelevel(cc,periodic,o_interior,overlap);
1145  init();
1146  }
1147 
1148  private:
1149 
1164  YaspGrid (std::array<std::vector<ctype>, dim> coords,
1165  std::bitset<dim> periodic,
1166  int overlap,
1168  std::array<int,dim> coarseSize,
1169  const YLoadBalance<dim>* lb = defaultLoadbalancer())
1170  : ccobj(comm), _torus(comm,tag,coarseSize,lb), leafIndexSet_(*this),
1171  _periodic(periodic), _coarseSize(coarseSize), _overlap(overlap),
1172  keep_ovlp(true), adaptRefCount(0), adaptActive(false)
1173  {
1174  // check whether YaspGrid has been given the correct template parameter
1175  static_assert(std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >::value,
1176  "YaspGrid coordinate container template parameter and given constructor values do not match!");
1177 
1178  if (!Dune::Yasp::checkIfMonotonous(coords))
1179  DUNE_THROW(Dune::GridError,"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1180 
1181  for (int i=0; i<dim; i++)
1182  _L[i] = coords[i][coords[i].size() - 1] - coords[i][0];
1183 
1184  _levels.resize(1);
1185 
1186  std::array<int,dim> o;
1187  std::fill(o.begin(), o.end(), 0);
1188  std::array<int,dim> o_interior(o);
1189  std::array<int,dim> s_interior(coarseSize);
1190 
1191  _torus.partition(_torus.rank(),o,coarseSize,o_interior,s_interior);
1192 
1193  // get offset by modifying o_interior according to overlap
1194  std::array<int,dim> offset(o_interior);
1195  for (int i=0; i<dim; i++)
1196  if ((periodic[i]) || (o_interior[i] > 0))
1197  offset[i] -= overlap;
1198 
1199  TensorProductCoordinates<ctype,dim> cc(coords, offset);
1200 
1201  // add level
1202  makelevel(cc,periodic,o_interior,overlap);
1203 
1204  init();
1205  }
1206 
1207  // the backup restore facility needs to be able to use above constructor
1208  friend struct BackupRestoreFacility<YaspGrid<dim,Coordinates> >;
1209 
1210  // do not copy this class
1211  YaspGrid(const YaspGrid&);
1212 
1213  public:
1214 
1218  int maxLevel() const
1219  {
1220  return _levels.size()-1;
1221  }
1222 
1224  void globalRefine (int refCount)
1225  {
1226  if (refCount < -maxLevel())
1227  DUNE_THROW(GridError, "Only " << maxLevel() << " levels left. " <<
1228  "Coarsening " << -refCount << " levels requested!");
1229 
1230  // If refCount is negative then coarsen the grid
1231  for (int k=refCount; k<0; k++)
1232  {
1233  // create an empty grid level
1234  YGridLevel empty;
1235  _levels.back() = empty;
1236  // reduce maxlevel
1237  _levels.pop_back();
1238 
1239  indexsets.pop_back();
1240  }
1241 
1242  // If refCount is positive refine the grid
1243  for (int k=0; k<refCount; k++)
1244  {
1245  // access to coarser grid level
1246  YGridLevel& cg = _levels[maxLevel()];
1247 
1248  std::bitset<dim> ovlp_low(0ULL), ovlp_up(0ULL);
1249  for (int i=0; i<dim; i++)
1250  {
1251  if (cg.overlap[0].dataBegin()->origin(i) > 0 || _periodic[i])
1252  ovlp_low[i] = true;
1253  if (cg.overlap[0].dataBegin()->max(i) + 1 < globalSize(i) || _periodic[i])
1254  ovlp_up[i] = true;
1255  }
1256 
1257  Coordinates newcont(cg.coords.refine(ovlp_low, ovlp_up, cg.overlapSize, keep_ovlp));
1258 
1259  int overlap = (keep_ovlp) ? 2*cg.overlapSize : cg.overlapSize;
1260 
1261  //determine new origin
1262  iTupel o_interior;
1263  for (int i=0; i<dim; i++)
1264  o_interior[i] = 2*cg.interior[0].dataBegin()->origin(i);
1265 
1266  // add level
1267  _levels.resize(_levels.size() + 1);
1268  makelevel(newcont,_periodic,o_interior,overlap);
1269 
1270  indexsets.push_back( std::make_shared<YaspIndexSet<const YaspGrid<dim,Coordinates>, false > >(*this,maxLevel()) );
1271  }
1272  }
1273 
1278  void refineOptions (bool keepPhysicalOverlap)
1279  {
1280  keep_ovlp = keepPhysicalOverlap;
1281  }
1282 
1294  bool mark( int refCount, const typename Traits::template Codim<0>::Entity & e )
1295  {
1296  assert(adaptActive == false);
1297  if (e.level() != maxLevel()) return false;
1298  adaptRefCount = std::max(adaptRefCount, refCount);
1299  return true;
1300  }
1301 
1308  int getMark ( const typename Traits::template Codim<0>::Entity &e ) const
1309  {
1310  return ( e.level() == maxLevel() ) ? adaptRefCount : 0;
1311  }
1312 
1314  bool adapt ()
1315  {
1316  globalRefine(adaptRefCount);
1317  return (adaptRefCount > 0);
1318  }
1319 
1321  bool preAdapt ()
1322  {
1323  adaptActive = true;
1324  adaptRefCount = comm().max(adaptRefCount);
1325  return (adaptRefCount < 0);
1326  }
1327 
1329  void postAdapt()
1330  {
1331  adaptActive = false;
1332  adaptRefCount = 0;
1333  }
1334 
1336  template<int cd, PartitionIteratorType pitype>
1337  typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lbegin (int level) const
1338  {
1339  return levelbegin<cd,pitype>(level);
1340  }
1341 
1343  template<int cd, PartitionIteratorType pitype>
1344  typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator lend (int level) const
1345  {
1346  return levelend<cd,pitype>(level);
1347  }
1348 
1350  template<int cd>
1351  typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lbegin (int level) const
1352  {
1353  return levelbegin<cd,All_Partition>(level);
1354  }
1355 
1357  template<int cd>
1358  typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator lend (int level) const
1359  {
1360  return levelend<cd,All_Partition>(level);
1361  }
1362 
1364  template<int cd, PartitionIteratorType pitype>
1365  typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafbegin () const
1366  {
1367  return levelbegin<cd,pitype>(maxLevel());
1368  }
1369 
1371  template<int cd, PartitionIteratorType pitype>
1372  typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator leafend () const
1373  {
1374  return levelend<cd,pitype>(maxLevel());
1375  }
1376 
1378  template<int cd>
1379  typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafbegin () const
1380  {
1381  return levelbegin<cd,All_Partition>(maxLevel());
1382  }
1383 
1385  template<int cd>
1386  typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator leafend () const
1387  {
1388  return levelend<cd,All_Partition>(maxLevel());
1389  }
1390 
1391  // \brief obtain Entity from EntitySeed. */
1392  template <typename Seed>
1393  typename Traits::template Codim<Seed::codimension>::Entity
1394  entity(const Seed& seed) const
1395  {
1396  const int codim = Seed::codimension;
1397  YGridLevelIterator g = begin(seed.impl().level());
1398 
1399  typedef typename Traits::template Codim<Seed::codimension>::Entity Entity;
1400  typedef YaspEntity<codim,dim,const YaspGrid> EntityImp;
1401  typedef typename YGrid::Iterator YIterator;
1402 
1403  return Entity(EntityImp(g,YIterator(g->overlapfront[codim],seed.impl().coord(),seed.impl().offset())));
1404  }
1405 
1407  int overlapSize (int level, int codim) const
1408  {
1409  YGridLevelIterator g = begin(level);
1410  return g->overlapSize;
1411  }
1412 
1414  int overlapSize (int codim) const
1415  {
1417  return g->overlapSize;
1418  }
1419 
1421  int ghostSize (int level, int codim) const
1422  {
1423  return 0;
1424  }
1425 
1427  int ghostSize (int codim) const
1428  {
1429  return 0;
1430  }
1431 
1433  int size (int level, int codim) const
1434  {
1435  YGridLevelIterator g = begin(level);
1436 
1437  // sum over all components of the codimension
1438  int count = 0;
1439  typedef typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator DAI;
1440  for (DAI it = g->overlapfront[codim].dataBegin(); it != g->overlapfront[codim].dataEnd(); ++it)
1441  count += it->totalsize();
1442 
1443  return count;
1444  }
1445 
1447  int size (int codim) const
1448  {
1449  return size(maxLevel(),codim);
1450  }
1451 
1453  int size (int level, GeometryType type) const
1454  {
1455  return (type.isCube()) ? size(level,dim-type.dim()) : 0;
1456  }
1457 
1459  int size (GeometryType type) const
1460  {
1461  return size(maxLevel(),type);
1462  }
1463 
1465  size_t numBoundarySegments () const
1466  {
1467  return nBSegments;
1468  }
1469 
1471  const Dune::FieldVector<ctype, dim>& domainSize () const {
1472  return _L;
1473  }
1474 
1479  template<class DataHandleImp, class DataType>
1481  {
1482  YaspCommunicateMeta<dim,dim>::comm(*this,data,iftype,dir,level);
1483  }
1484 
1489  template<class DataHandleImp, class DataType>
1491  {
1492  YaspCommunicateMeta<dim,dim>::comm(*this,data,iftype,dir,this->maxLevel());
1493  }
1494 
1499  template<class DataHandle, int codim>
1500  void communicateCodim (DataHandle& data, InterfaceType iftype, CommunicationDirection dir, int level) const
1501  {
1502  // check input
1503  if (!data.contains(dim,codim)) return; // should have been checked outside
1504 
1505  // data types
1506  typedef typename DataHandle::DataType DataType;
1507 
1508  // access to grid level
1509  YGridLevelIterator g = begin(level);
1510 
1511  // find send/recv lists or throw error
1512  const YGridList<Coordinates>* sendlist = 0;
1513  const YGridList<Coordinates>* recvlist = 0;
1514 
1516  {
1517  sendlist = &g->send_interiorborder_interiorborder[codim];
1518  recvlist = &g->recv_interiorborder_interiorborder[codim];
1519  }
1520  if (iftype==InteriorBorder_All_Interface)
1521  {
1522  sendlist = &g->send_interiorborder_overlapfront[codim];
1523  recvlist = &g->recv_overlapfront_interiorborder[codim];
1524  }
1526  {
1527  sendlist = &g->send_overlap_overlapfront[codim];
1528  recvlist = &g->recv_overlapfront_overlap[codim];
1529  }
1530  if (iftype==All_All_Interface)
1531  {
1532  sendlist = &g->send_overlapfront_overlapfront[codim];
1533  recvlist = &g->recv_overlapfront_overlapfront[codim];
1534  }
1535 
1536  // change communication direction?
1537  if (dir==BackwardCommunication)
1538  std::swap(sendlist,recvlist);
1539 
1540  int cnt;
1541 
1542  // Size computation (requires communication if variable size)
1543  std::vector<int> send_size(sendlist->size(),-1); // map rank to total number of objects (of type DataType) to be sent
1544  std::vector<int> recv_size(recvlist->size(),-1); // map rank to total number of objects (of type DataType) to be recvd
1545  std::vector<size_t*> send_sizes(sendlist->size(),static_cast<size_t*>(0)); // map rank to array giving number of objects per entity to be sent
1546  std::vector<size_t*> recv_sizes(recvlist->size(),static_cast<size_t*>(0)); // map rank to array giving number of objects per entity to be recvd
1547 
1548  // define type to iterate over send and recv lists
1549  typedef typename YGridList<Coordinates>::Iterator ListIt;
1550 
1551  if (data.fixedSize(dim,codim))
1552  {
1553  // fixed size: just take a dummy entity, size can be computed without communication
1554  cnt=0;
1555  for (ListIt is=sendlist->begin(); is!=sendlist->end(); ++is)
1556  {
1557  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1559  send_size[cnt] = is->grid.totalsize() * data.size(*it);
1560  cnt++;
1561  }
1562  cnt=0;
1563  for (ListIt is=recvlist->begin(); is!=recvlist->end(); ++is)
1564  {
1565  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1567  recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1568  cnt++;
1569  }
1570  }
1571  else
1572  {
1573  // variable size case: sender side determines the size
1574  cnt=0;
1575  for (ListIt is=sendlist->begin(); is!=sendlist->end(); ++is)
1576  {
1577  // allocate send buffer for sizes per entitiy
1578  size_t *buf = new size_t[is->grid.totalsize()];
1579  send_sizes[cnt] = buf;
1580 
1581  // loop over entities and ask for size
1582  int i=0; size_t n=0;
1583  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1585  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1586  itend(YaspLevelIterator<codim,All_Partition,GridImp>(g, typename YGrid::Iterator(is->yg,true)));
1587  for ( ; it!=itend; ++it)
1588  {
1589  buf[i] = data.size(*it);
1590  n += buf[i];
1591  i++;
1592  }
1593 
1594  // now we know the size for this rank
1595  send_size[cnt] = n;
1596 
1597  // hand over send request to torus class
1598  torus().send(is->rank,buf,is->grid.totalsize()*sizeof(size_t));
1599  cnt++;
1600  }
1601 
1602  // allocate recv buffers for sizes and store receive request
1603  cnt=0;
1604  for (ListIt is=recvlist->begin(); is!=recvlist->end(); ++is)
1605  {
1606  // allocate recv buffer
1607  size_t *buf = new size_t[is->grid.totalsize()];
1608  recv_sizes[cnt] = buf;
1609 
1610  // hand over recv request to torus class
1611  torus().recv(is->rank,buf,is->grid.totalsize()*sizeof(size_t));
1612  cnt++;
1613  }
1614 
1615  // exchange all size buffers now
1616  torus().exchange();
1617 
1618  // release send size buffers
1619  cnt=0;
1620  for (ListIt is=sendlist->begin(); is!=sendlist->end(); ++is)
1621  {
1622  delete[] send_sizes[cnt];
1623  send_sizes[cnt] = 0;
1624  cnt++;
1625  }
1626 
1627  // process receive size buffers
1628  cnt=0;
1629  for (ListIt is=recvlist->begin(); is!=recvlist->end(); ++is)
1630  {
1631  // get recv buffer
1632  size_t *buf = recv_sizes[cnt];
1633 
1634  // compute total size
1635  size_t n=0;
1636  for (int i=0; i<is->grid.totalsize(); ++i)
1637  n += buf[i];
1638 
1639  // ... and store it
1640  recv_size[cnt] = n;
1641  ++cnt;
1642  }
1643  }
1644 
1645 
1646  // allocate & fill the send buffers & store send request
1647  std::vector<DataType*> sends(sendlist->size(), static_cast<DataType*>(0)); // store pointers to send buffers
1648  cnt=0;
1649  for (ListIt is=sendlist->begin(); is!=sendlist->end(); ++is)
1650  {
1651  // allocate send buffer
1652  DataType *buf = new DataType[send_size[cnt]];
1653 
1654  // remember send buffer
1655  sends[cnt] = buf;
1656 
1657  // make a message buffer
1658  MessageBuffer<DataType> mb(buf);
1659 
1660  // fill send buffer; iterate over cells in intersection
1661  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1663  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1664  itend(YaspLevelIterator<codim,All_Partition,GridImp>(g, typename YGrid::Iterator(is->yg,true)));
1665  for ( ; it!=itend; ++it)
1666  data.gather(mb,*it);
1667 
1668  // hand over send request to torus class
1669  torus().send(is->rank,buf,send_size[cnt]*sizeof(DataType));
1670  cnt++;
1671  }
1672 
1673  // allocate recv buffers and store receive request
1674  std::vector<DataType*> recvs(recvlist->size(),static_cast<DataType*>(0)); // store pointers to send buffers
1675  cnt=0;
1676  for (ListIt is=recvlist->begin(); is!=recvlist->end(); ++is)
1677  {
1678  // allocate recv buffer
1679  DataType *buf = new DataType[recv_size[cnt]];
1680 
1681  // remember recv buffer
1682  recvs[cnt] = buf;
1683 
1684  // hand over recv request to torus class
1685  torus().recv(is->rank,buf,recv_size[cnt]*sizeof(DataType));
1686  cnt++;
1687  }
1688 
1689  // exchange all buffers now
1690  torus().exchange();
1691 
1692  // release send buffers
1693  cnt=0;
1694  for (ListIt is=sendlist->begin(); is!=sendlist->end(); ++is)
1695  {
1696  delete[] sends[cnt];
1697  sends[cnt] = 0;
1698  cnt++;
1699  }
1700 
1701  // process receive buffers and delete them
1702  cnt=0;
1703  for (ListIt is=recvlist->begin(); is!=recvlist->end(); ++is)
1704  {
1705  // get recv buffer
1706  DataType *buf = recvs[cnt];
1707 
1708  // make a message buffer
1709  MessageBuffer<DataType> mb(buf);
1710 
1711  // copy data from receive buffer; iterate over cells in intersection
1712  if (data.fixedSize(dim,codim))
1713  {
1714  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1716  size_t n=data.size(*it);
1717  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1718  itend(YaspLevelIterator<codim,All_Partition,GridImp>(g, typename YGrid::Iterator(is->yg,true)));
1719  for ( ; it!=itend; ++it)
1720  data.scatter(mb,*it,n);
1721  }
1722  else
1723  {
1724  int i=0;
1725  size_t *sbuf = recv_sizes[cnt];
1726  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1728  typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1729  itend(YaspLevelIterator<codim,All_Partition,GridImp>(g, typename YGrid::Iterator(is->yg,true)));
1730  for ( ; it!=itend; ++it)
1731  data.scatter(mb,*it,sbuf[i++]);
1732  delete[] sbuf;
1733  }
1734 
1735  // delete buffer
1736  delete[] buf; // hier krachts !
1737  cnt++;
1738  }
1739  }
1740 
1741  // The new index sets from DDM 11.07.2005
1742  const typename Traits::GlobalIdSet& globalIdSet() const
1743  {
1744  return theglobalidset;
1745  }
1746 
1747  const typename Traits::LocalIdSet& localIdSet() const
1748  {
1749  return theglobalidset;
1750  }
1751 
1752  const typename Traits::LevelIndexSet& levelIndexSet(int level) const
1753  {
1754  if (level<0 || level>maxLevel()) DUNE_THROW(RangeError, "level out of range");
1755  return *(indexsets[level]);
1756  }
1757 
1758  const typename Traits::LeafIndexSet& leafIndexSet() const
1759  {
1760  return leafIndexSet_;
1761  }
1762 
1766  {
1767  return ccobj;
1768  }
1769 
1770  private:
1771 
1772  // number of boundary segments of the level 0 grid
1773  int nBSegments;
1774 
1775  // Index classes need access to the real entity
1776  friend class Dune::YaspIndexSet<const Dune::YaspGrid<dim, Coordinates>, true >;
1777  friend class Dune::YaspIndexSet<const Dune::YaspGrid<dim, Coordinates>, false >;
1778  friend class Dune::YaspGlobalIdSet<const Dune::YaspGrid<dim, Coordinates> >;
1779  friend class Dune::YaspPersistentContainerIndex<const Dune::YaspGrid<dim, Coordinates> >;
1780 
1781  friend class Dune::YaspIntersectionIterator<const Dune::YaspGrid<dim, Coordinates> >;
1782  friend class Dune::YaspIntersection<const Dune::YaspGrid<dim, Coordinates> >;
1783  friend class Dune::YaspEntity<0, dim, const Dune::YaspGrid<dim, Coordinates> >;
1784 
1785  template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1786  friend class Entity;
1787 
1788  template<class DT>
1789  class MessageBuffer {
1790  public:
1791  // Constructor
1792  MessageBuffer (DT *p)
1793  {
1794  a=p;
1795  i=0;
1796  j=0;
1797  }
1798 
1799  // write data to message buffer, acts like a stream !
1800  template<class Y>
1801  void write (const Y& data)
1802  {
1803  static_assert(( std::is_same<DT,Y>::value ), "DataType mismatch");
1804  a[i++] = data;
1805  }
1806 
1807  // read data from message buffer, acts like a stream !
1808  template<class Y>
1809  void read (Y& data) const
1810  {
1811  static_assert(( std::is_same<DT,Y>::value ), "DataType mismatch");
1812  data = a[j++];
1813  }
1814 
1815  private:
1816  DT *a;
1817  int i;
1818  mutable int j;
1819  };
1820 
1822  template<int cd, PartitionIteratorType pitype>
1823  YaspLevelIterator<cd,pitype,GridImp> levelbegin (int level) const
1824  {
1825  YGridLevelIterator g = begin(level);
1826  if (level<0 || level>maxLevel()) DUNE_THROW(RangeError, "level out of range");
1827 
1828  if (pitype==Interior_Partition)
1829  return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].begin());
1830  if (pitype==InteriorBorder_Partition)
1831  return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].begin());
1832  if (pitype==Overlap_Partition)
1833  return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].begin());
1834  if (pitype<=All_Partition)
1835  return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].begin());
1836  if (pitype==Ghost_Partition)
1837  return levelend <cd, pitype> (level);
1838 
1839  DUNE_THROW(GridError, "YaspLevelIterator with this codim or partition type not implemented");
1840  }
1841 
1843  template<int cd, PartitionIteratorType pitype>
1844  YaspLevelIterator<cd,pitype,GridImp> levelend (int level) const
1845  {
1846  YGridLevelIterator g = begin(level);
1847  if (level<0 || level>maxLevel()) DUNE_THROW(RangeError, "level out of range");
1848 
1849  if (pitype==Interior_Partition)
1850  return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].end());
1851  if (pitype==InteriorBorder_Partition)
1852  return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].end());
1853  if (pitype==Overlap_Partition)
1854  return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].end());
1855  if (pitype<=All_Partition || pitype == Ghost_Partition)
1856  return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].end());
1857 
1858  DUNE_THROW(GridError, "YaspLevelIterator with this codim or partition type not implemented");
1859  }
1860 
1862 
1863  Torus<CollectiveCommunicationType,dim> _torus;
1864 
1865  std::vector< std::shared_ptr< YaspIndexSet<const YaspGrid<dim,Coordinates>, false > > > indexsets;
1866  YaspIndexSet<const YaspGrid<dim,Coordinates>, true> leafIndexSet_;
1867  YaspGlobalIdSet<const YaspGrid<dim,Coordinates> > theglobalidset;
1868 
1869  Dune::FieldVector<ctype, dim> _L;
1870  iTupel _s;
1871  std::bitset<dim> _periodic;
1872  iTupel _coarseSize;
1873  ReservedVector<YGridLevel,32> _levels;
1874  int _overlap;
1875  bool keep_ovlp;
1876  int adaptRefCount;
1877  bool adaptActive;
1878  };
1879 
1881 
1882  template <int d, class CC>
1883  std::ostream& operator<< (std::ostream& s, const YaspGrid<d,CC>& grid)
1884  {
1885  int rank = grid.torus().rank();
1886 
1887  s << "[" << rank << "]:" << " YaspGrid maxlevel=" << grid.maxLevel() << std::endl;
1888 
1889  s << "Printing the torus: " <<std::endl;
1890  s << grid.torus() << std::endl;
1891 
1892  for (typename YaspGrid<d,CC>::YGridLevelIterator g=grid.begin(); g!=grid.end(); ++g)
1893  {
1894  s << "[" << rank << "]: " << std::endl;
1895  s << "[" << rank << "]: " << "==========================================" << std::endl;
1896  s << "[" << rank << "]: " << "level=" << g->level() << std::endl;
1897 
1898  for (int codim = 0; codim < d + 1; ++codim)
1899  {
1900  s << "[" << rank << "]: " << "overlapfront[" << codim << "]: " << g->overlapfront[codim] << std::endl;
1901  s << "[" << rank << "]: " << "overlap[" << codim << "]: " << g->overlap[codim] << std::endl;
1902  s << "[" << rank << "]: " << "interiorborder[" << codim << "]: " << g->interiorborder[codim] << std::endl;
1903  s << "[" << rank << "]: " << "interior[" << codim << "]: " << g->interior[codim] << std::endl;
1904 
1905  typedef typename YGridList<CC>::Iterator I;
1906  for (I i=g->send_overlapfront_overlapfront[codim].begin();
1907  i!=g->send_overlapfront_overlapfront[codim].end(); ++i)
1908  s << "[" << rank << "]: " << " s_of_of[" << codim << "] to rank "
1909  << i->rank << " " << i->grid << std::endl;
1910 
1911  for (I i=g->recv_overlapfront_overlapfront[codim].begin();
1912  i!=g->recv_overlapfront_overlapfront[codim].end(); ++i)
1913  s << "[" << rank << "]: " << " r_of_of[" << codim << "] to rank "
1914  << i->rank << " " << i->grid << std::endl;
1915 
1916  for (I i=g->send_overlap_overlapfront[codim].begin();
1917  i!=g->send_overlap_overlapfront[codim].end(); ++i)
1918  s << "[" << rank << "]: " << " s_o_of[" << codim << "] to rank "
1919  << i->rank << " " << i->grid << std::endl;
1920 
1921  for (I i=g->recv_overlapfront_overlap[codim].begin();
1922  i!=g->recv_overlapfront_overlap[codim].end(); ++i)
1923  s << "[" << rank << "]: " << " r_of_o[" << codim << "] to rank "
1924  << i->rank << " " << i->grid << std::endl;
1925 
1926  for (I i=g->send_interiorborder_interiorborder[codim].begin();
1927  i!=g->send_interiorborder_interiorborder[codim].end(); ++i)
1928  s << "[" << rank << "]: " << " s_ib_ib[" << codim << "] to rank "
1929  << i->rank << " " << i->grid << std::endl;
1930 
1931  for (I i=g->recv_interiorborder_interiorborder[codim].begin();
1932  i!=g->recv_interiorborder_interiorborder[codim].end(); ++i)
1933  s << "[" << rank << "]: " << " r_ib_ib[" << codim << "] to rank "
1934  << i->rank << " " << i->grid << std::endl;
1935 
1936  for (I i=g->send_interiorborder_overlapfront[codim].begin();
1937  i!=g->send_interiorborder_overlapfront[codim].end(); ++i)
1938  s << "[" << rank << "]: " << " s_ib_of[" << codim << "] to rank "
1939  << i->rank << " " << i->grid << std::endl;
1940 
1941  for (I i=g->recv_overlapfront_interiorborder[codim].begin();
1942  i!=g->recv_overlapfront_interiorborder[codim].end(); ++i)
1943  s << "[" << rank << "]: " << " r_of_ib[" << codim << "] to rank "
1944  << i->rank << " " << i->grid << std::endl;
1945  }
1946  }
1947 
1948  s << std::endl;
1949 
1950  return s;
1951  }
1952 
1953  namespace Capabilities
1954  {
1955 
1963  template<int dim, class Coordinates>
1964  struct hasBackupRestoreFacilities< YaspGrid<dim, Coordinates> >
1965  {
1966  static const bool v = true;
1967  };
1968 
1972  template<int dim, class Coordinates>
1973  struct hasSingleGeometryType< YaspGrid<dim, Coordinates> >
1974  {
1975  static const bool v = true;
1976  static const unsigned int topologyId = GeometryTypes::cube(dim).id();
1977  };
1978 
1982  template<int dim, class Coordinates>
1983  struct isCartesian< YaspGrid<dim, Coordinates> >
1984  {
1985  static const bool v = true;
1986  };
1987 
1991  template<int dim, class Coordinates, int codim>
1992  struct hasEntity< YaspGrid<dim, Coordinates>, codim>
1993  {
1994  static const bool v = true;
1995  };
1996 
2001  template<int dim, class Coordinates, int codim>
2002  struct hasEntityIterator<YaspGrid<dim, Coordinates>, codim>
2003  {
2004  static const bool v = true;
2005  };
2006 
2010  template<int dim, int codim, class Coordinates>
2011  struct canCommunicate< YaspGrid< dim, Coordinates>, codim >
2012  {
2013  static const bool v = true;
2014  };
2015 
2019  template<int dim, class Coordinates>
2020  struct isLevelwiseConforming< YaspGrid<dim, Coordinates> >
2021  {
2022  static const bool v = true;
2023  };
2024 
2028  template<int dim, class Coordinates>
2029  struct isLeafwiseConforming< YaspGrid<dim, Coordinates> >
2030  {
2031  static const bool v = true;
2032  };
2033 
2034  }
2035 
2036 } // end namespace
2037 
2038 // Include the specialization of the StructuredGridFactory class for YaspGrid
2040 // Include the specialization of the BackupRestoreFacility class for YaspGrid
2042 
2043 #endif
Dune::Capabilities::canCommunicate::v
static const bool v
Definition: common/capabilities.hh:96
Dune::YaspGrid::YaspGrid
YaspGrid(std::array< std::vector< ctype >, dim > coords, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Standard constructor for a tensorproduct YaspGrid.
Definition: yaspgrid.hh:1047
Dune::YGridComponent
Definition: ygrid.hh:72
yaspgridpersistentcontainer.hh
Specialization of the PersistentContainer for YaspGrid.
Dune::YaspGridFamily::CCType
CollectiveCommunication< No_Comm > CCType
Definition: yaspgrid.hh:94
Dune::Partitions::overlap
constexpr Overlap overlap
PartitionSet for the overlap partition.
Definition: partitionset.hh:276
Dune::Torus::dims
const iTupel & dims() const
return dimensions of torus
Definition: torus.hh:110
Dune::Torus::sendend
ProcListIterator sendend() const
end of send list
Definition: torus.hh:341
Dune::YaspGrid
[ provides Dune::Grid ]
Definition: yaspgrid.hh:58
uint8_t
unsigned char uint8_t
Definition: yaspgrid.hh:16
Dune::YGridList
implements a collection of multiple std::deque<Intersection> Intersections with neighboring processor...
Definition: ygrid.hh:820
Dune::EquidistantCoordinates
Container for equidistant coordinates in a YaspGrid.
Definition: coordinates.hh:26
Dune::Capabilities::isLevelwiseConforming::v
static const bool v
Definition: common/capabilities.hh:105
Dune::YaspGrid::intersections
void intersections(const YGridComponent< Coordinates > &sendgrid, const YGridComponent< Coordinates > &recvgrid, std::deque< Intersection > &sendlist, std::deque< Intersection > &recvlist)
Construct list of intersections with neighboring processors.
Definition: yaspgrid.hh:563
Dune::Capabilities::hasEntityIterator
specialize with 'true' for all codims that a grid provides an iterator for (default=false)
Definition: common/capabilities.hh:71
Dune::YaspGrid::size
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: yaspgrid.hh:1459
Dune::YaspGrid::leafend
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1372
Dune::YaspGrid::Traits
YaspGridFamily< dim, Coordinates >::Traits Traits
Definition: yaspgrid.hh:719
Dune::YaspGrid::YGridLevelIterator
ReservedVector< YGridLevel, 32 >::const_iterator YGridLevelIterator
Iterator over the grid levels.
Definition: yaspgrid.hh:293
Dune::InteriorBorder_Partition
@ InteriorBorder_Partition
interior and border entities
Definition: gridenums.hh:136
indexidset.hh
Provides base classes for index and id sets.
Dune::YaspGrid::init
void init()
Definition: yaspgrid.hh:679
Dune::YaspGrid::makelevel
void makelevel(const Coordinates &coords, std::bitset< dim > periodic, iTupel o_interior, int overlap)
Make a new YGridLevel structure.
Definition: yaspgrid.hh:330
Dune::Capabilities::isLeafwiseConforming
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:112
Dune::YaspGrid::begin
YGridLevelIterator begin(int i) const
return iterator pointing to given level
Definition: yaspgrid.hh:302
Dune::YGrid
implements a collection of YGridComponents which form a codimension Entities of given codimension c n...
Definition: ygrid.hh:548
Dune::Capabilities::hasEntity::v
static const bool v
Definition: common/capabilities.hh:57
Dune::YaspGrid::levelSize
iTupel levelSize(int l) const
return size vector of the grid (in cells) on level l
Definition: yaspgrid.hh:273
Dune::YaspGrid::getRefineOption
bool getRefineOption() const
Definition: yaspgrid.hh:287
yaspgridentityseed.hh
The YaspEntitySeed class.
Dune::Torus< CollectiveCommunicationType, dim >
Dune::YaspGlobalIdSet
persistent, globally unique Ids
Definition: yaspgrid.hh:67
Dune::Overlap_Partition
@ Overlap_Partition
interior, border, and overlap entities
Definition: gridenums.hh:137
Dune::Yasp::sizeArray
std::array< int, d > sizeArray(const std::array< std::vector< ct >, d > &v)
Definition: ygrid.hh:27
Dune::YaspGrid::postAdapt
void postAdapt()
clean up some markers
Definition: yaspgrid.hh:1329
Dune::Capabilities::canCommunicate
specialize with 'true' for all codims that a grid can communicate data on (default=false)
Definition: common/capabilities.hh:94
Dune::YaspGrid::PersistentIndexType
bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim > PersistentIndexType
Definition: yaspgrid.hh:714
Dune::Torus::send
void send(int rank, void *buffer, int size) const
store a send request; buffers are sent in order; handles also local requests with memcpy
Definition: torus.hh:359
Dune::Torus::coord
iTupel coord() const
return own coordinates
Definition: torus.hh:98
Dune::Capabilities::hasEntity
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:55
Dune::Capabilities::hasBackupRestoreFacilities::v
static const bool v
Definition: common/capabilities.hh:123
Dune::GridError
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
Dune::operator<<
std::ostream & operator<<(std::ostream &out, const PartitionType &type)
write a PartitionType to a stream
Definition: gridenums.hh:70
Dune::InteriorBorder_InteriorBorder_Interface
@ InteriorBorder_InteriorBorder_Interface
send/receive interior and border entities
Definition: gridenums.hh:85
Dune::YaspHierarchicIterator
YaspHierarchicIterator enables iteration over son entities of codim 0.
Definition: yaspgrid.hh:65
Dune::YaspGrid::leafbegin
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1379
Dune::YLoadBalanceDefault
Implement the default load balance strategy of yaspgrid.
Definition: partitioning.hh:34
Dune::YaspLevelIterator
Iterates over entities of one grid level.
Definition: yaspgrid.hh:62
Dune::YaspGrid::communicate
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1480
Dune::Capabilities::hasSingleGeometryType
Specialize with 'true' for if the codimension 0 entity of the grid has only one possible geometry typ...
Definition: common/capabilities.hh:24
Dune::Capabilities::isLeafwiseConforming::v
static const bool v
Definition: common/capabilities.hh:114
Dune::DefaultLevelGridViewTraits
Definition: defaultgridview.hh:23
Dune::IdSet
Id Set Interface.
Definition: common/grid.hh:347
Dune::YaspGrid::globalIdSet
const Traits::GlobalIdSet & globalIdSet() const
Definition: yaspgrid.hh:1742
Dune::BackwardCommunication
@ BackwardCommunication
reverse communication direction
Definition: gridenums.hh:170
Dune::YaspGrid::mark
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Marks an entity to be refined/coarsened in a subsequent adapt.
Definition: yaspgrid.hh:1294
Dune::YaspGrid::globalSize
iTupel globalSize() const
return number of cells on finest level on all processors
Definition: yaspgrid.hh:261
Dune::CommDataHandleIF
CommDataHandleIF describes the features of a data handle for communication in parallel runs using the...
Definition: datahandleif.hh:75
Dune::Intersection
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/grid.hh:344
Dune::yaspgrid_dim_bits
const int yaspgrid_dim_bits
Definition: yaspgrid.hh:51
Dune::YaspGrid::leafIndexSet
const Traits::LeafIndexSet & leafIndexSet() const
Definition: yaspgrid.hh:1758
Dune::YaspGrid::YaspGrid
YaspGrid(Dune::FieldVector< ctype, dim > L, std::array< int, dim > s, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:904
Dune::Torus::exchange
void exchange() const
exchange messages stored in request buffers; clear request buffers afterwards
Definition: torus.hh:385
Dune::YaspGrid::boundarysegmentssize
void boundarysegmentssize()
Definition: yaspgrid.hh:685
Dune::YaspGrid::GridImp
const typedef YaspGrid< dim, Coordinates > GridImp
Definition: yaspgrid.hh:677
Dune::YaspGrid::CollectiveCommunicationType
CollectiveCommunication< No_Comm > CollectiveCommunicationType
Definition: yaspgrid.hh:183
Dune::YGridList::Iterator
Definition: ygrid.hh:843
Dune::TensorProductCoordinates
Coordinate container for a tensor product YaspGrid.
Definition: coordinates.hh:242
Dune::YaspEntity
Definition: yaspgrid.hh:60
Dune::Interior_Partition
@ Interior_Partition
only interior entities
Definition: gridenums.hh:135
Dune::Torus::partition
double partition(int rank, iTupel origin_in, iTupel size_in, iTupel &origin_out, iTupel &size_out) const
partition the given grid onto the torus and return the piece of the process with given rank; returns ...
Definition: torus.hh:237
Dune::YaspGridFamily::Traits
GridTraits< dim, dim, Dune::YaspGrid< dim, Coordinates >, YaspGeometry, YaspEntity, YaspLevelIterator, YaspIntersection, YaspIntersection, YaspIntersectionIterator, YaspIntersectionIterator, YaspHierarchicIterator, YaspLevelIterator, YaspIndexSet< const YaspGrid< dim, Coordinates >, false >, YaspIndexSet< const YaspGrid< dim, Coordinates >, true >, YaspGlobalIdSet< const YaspGrid< dim, Coordinates > >, bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim >, YaspGlobalIdSet< const YaspGrid< dim, Coordinates > >, bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim >, CCType, DefaultLevelGridViewTraits, DefaultLeafGridViewTraits, YaspEntitySeed > Traits
Definition: yaspgrid.hh:117
Dune::Torus::recvend
ProcListIterator recvend() const
last process in receive list
Definition: torus.hh:353
Dune::YaspGrid::defaultLoadbalancer
static const YLoadBalanceDefault< dim > * defaultLoadbalancer()
Definition: yaspgrid.hh:316
Dune::All_All_Interface
@ All_All_Interface
send all and receive all entities
Definition: gridenums.hh:89
Dune::All_Partition
@ All_Partition
all entities
Definition: gridenums.hh:139
Dune::Torus::sendbegin
ProcListIterator sendbegin() const
first process in send list
Definition: torus.hh:335
Dune::YaspGrid::LeafIndexSetType
YaspIndexSet< YaspGrid< dim, Coordinates >, true > LeafIndexSetType
Definition: yaspgrid.hh:723
Dune::YaspGrid::globalRefine
void globalRefine(int refCount)
refine the grid refCount times.
Definition: yaspgrid.hh:1224
Dune::Torus::recvbegin
ProcListIterator recvbegin() const
first process in receive list
Definition: torus.hh:347
Dune::YaspGrid::size
int size(int codim) const
number of leaf entities per codim in this process
Definition: yaspgrid.hh:1447
Dune::YaspGrid::localIdSet
const Traits::LocalIdSet & localIdSet() const
Definition: yaspgrid.hh:1747
Dune::YaspGrid::lend
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lend(int level) const
Iterator to one past the last entity of given codim on level for partition type.
Definition: yaspgrid.hh:1344
Dune::YaspGridFamily
Definition: yaspgrid.hh:89
coordinates.hh
This provides container classes for the coordinates to be used in YaspGrid Upon implementation of the...
Dune::Yasp::checkIfMonotonous
bool checkIfMonotonous(const std::array< std::vector< ctype >, dim > &coords)
Definition: coordinates.hh:370
yaspgrididset.hh
Dune::Capabilities::hasEntityIterator::v
static const bool v
Definition: common/capabilities.hh:73
Dune::EquidistantOffsetCoordinates::origin
ct origin(int d) const
Definition: coordinates.hh:183
yaspgridentity.hh
the YaspEntity class and its specializations
Dune::YaspGrid::YaspGrid
YaspGrid(Dune::FieldVector< ctype, dim > lowerleft, Dune::FieldVector< ctype, dim > upperright, std::array< int, dim > s, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:976
Dune::YaspGrid::getMark
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
returns adaptation mark for given entity
Definition: yaspgrid.hh:1308
structuredyaspgridfactory.hh
Specialization of the StructuredGridFactory class for YaspGrid.
Dune::YaspGrid::comm
const CollectiveCommunicationType & comm() const
return a collective communication object
Definition: yaspgrid.hh:1765
Dune::YaspGrid::torus
const Torus< CollectiveCommunicationType, dim > & torus() const
return reference to torus
Definition: yaspgrid.hh:249
Dune::YaspIntersection
YaspIntersection provides data about intersection with neighboring codim 0 entities.
Definition: yaspgrid.hh:64
Dune::yaspgrid_level_bits
const int yaspgrid_level_bits
Definition: yaspgrid.hh:52
Dune::Ghost_Partition
@ Ghost_Partition
only ghost entities
Definition: gridenums.hh:140
Dune::YaspGrid::refineOptions
void refineOptions(bool keepPhysicalOverlap)
set options for refinement
Definition: yaspgrid.hh:1278
Dune::YaspPersistentContainerIndex
Definition: yaspgrid.hh:68
capabilities.hh
A set of traits classes to store static information about grid implementation.
Dune::YaspGrid::globalSize
int globalSize(int i) const
return number of cells on finest level in given direction on all processors
Definition: yaspgrid.hh:255
Dune::Torus::rank
int rank() const
return own rank
Definition: torus.hh:92
backuprestore.hh
Dune::YGridComponent::intersection
YGridComponent< Coordinates > intersection(const YGridComponent< Coordinates > &r) const
Return YGridComponent of supergrid of self which is the intersection of self and another YGridCompone...
Definition: ygrid.hh:269
Dune::YaspGrid::end
YGridLevelIterator end() const
return iterator pointing to one past the finest level
Definition: yaspgrid.hh:310
Dune::YaspGrid::isPeriodic
bool isPeriodic(int i) const
return whether the grid is periodic in direction i
Definition: yaspgrid.hh:282
Dune::Overlap_OverlapFront_Interface
@ Overlap_OverlapFront_Interface
send overlap, receive overlap and front entities
Definition: gridenums.hh:87
Dune::YGridList::size
int size() const
return the size of the container, this is the sum of the sizes of all deques
Definition: ygrid.hh:951
yaspgridintersectioniterator.hh
The YaspIntersectionIterator class.
Dune::BackupRestoreFacility
facility for writing and reading grids
Definition: common/backuprestore.hh:40
Dune::GridDefaultImplementation
Definition: common/geometry.hh:24
yaspgridgeometry.hh
The YaspGeometry class and its specializations.
Dune::YaspGrid::levelIndexSet
const Traits::LevelIndexSet & levelIndexSet(int level) const
Definition: yaspgrid.hh:1752
Dune::YGrid::Iterator
Iterator over a collection o YGrids A YGrid::Iterator is the heart of an entity in YaspGrid.
Definition: ygrid.hh:591
Dune::YaspIndexSet
Implementation of Level- and LeafIndexSets for YaspGrid.
Definition: yaspgrid.hh:66
Dune::YaspGrid::levelSize
int levelSize(int l, int i) const
return size of the grid (in cells) on level l in direction i
Definition: yaspgrid.hh:267
ygrid.hh
This provides a YGrid, the elemental component of the yaspgrid implementation.
Dune::YaspIntersectionIterator
YaspIntersectionIterator enables iteration over intersections with neighboring codim 0 entities.
Definition: yaspgrid.hh:63
Dune::YaspGrid::overlapSize
int overlapSize(int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1414
Dune::YaspGrid::leafend
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1386
Dune::YaspGrid::preAdapt
bool preAdapt()
returns true, if the grid will be coarsened
Definition: yaspgrid.hh:1321
Dune::Capabilities::isLevelwiseConforming
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:103
Dune::GridTraits
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1064
Dune::InteriorBorder_All_Interface
@ InteriorBorder_All_Interface
send interior and border, receive all entities
Definition: gridenums.hh:86
Dune::YaspGrid::YaspGrid
YaspGrid(const Coordinates &coordinates, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:733
backuprestore.hh
Dune::YLoadBalance
a base class for the yaspgrid partitioning strategy The name might be irritating. It will probably ch...
Definition: partitioning.hh:23
Dune::YaspEntitySeed
Describes the minimal information necessary to create a fully functional YaspEntity.
Definition: yaspgrid.hh:61
Dune::YaspGrid::lbegin
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lbegin(int level) const
one past the end on this level
Definition: yaspgrid.hh:1337
Dune::YaspGrid::size
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: yaspgrid.hh:1453
Dune::YGridList::begin
Iterator begin() const
return iterator pointing to the begin of the container
Definition: ygrid.hh:921
Dune::Partitions::front
constexpr Front front
PartitionSet for the front partition.
Definition: partitionset.hh:279
Dune::YGridComponent::move
YGridComponent< Coordinates > move(iTupel v) const
return grid moved by the vector v
Definition: ygrid.hh:261
yaspgridindexsets.hh
level-wise, non-persistent, consecutive indices for YaspGrid
Dune::DefaultLeafGridViewTraits
Definition: defaultgridview.hh:206
Dune::YaspGrid::numBoundarySegments
size_t numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: yaspgrid.hh:1465
Dune::YaspGrid::lbegin
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lbegin(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1351
datahandleif.hh
Describes the parallel communication interface class for MessageBuffers and DataHandles.
Dune::Partitions::interior
constexpr Interior interior
PartitionSet for the interior partition.
Definition: partitionset.hh:270
Dune::Capabilities::isCartesian
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: common/capabilities.hh:45
Dune::YaspGrid::ctype
Coordinates::ctype ctype
Type used for coordinates.
Definition: yaspgrid.hh:179
Dune::YaspGrid::begin
YGridLevelIterator begin() const
return iterator pointing to coarsest level
Definition: yaspgrid.hh:296
Dune::YaspGrid::adapt
bool adapt()
map adapt to global refine
Definition: yaspgrid.hh:1314
Dune::YaspGrid::maxLevel
int maxLevel() const
Definition: yaspgrid.hh:1218
yaspgridleveliterator.hh
The YaspLevelIterator class.
yaspgridhierarchiciterator.hh
Dune::Capabilities::hasSingleGeometryType::v
static const bool v
Definition: common/capabilities.hh:26
Dune::YaspGrid::ghostSize
int ghostSize(int level, int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1421
Dune::Torus::neighbors
int neighbors() const
return the number of neighbors, which is
Definition: torus.hh:201
torus.hh
This file provides the infrastructure for toroidal communication in YaspGrid.
Dune::YaspGrid::lend
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lend(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1358
Dune::YaspGrid::communicateCodim
void communicateCodim(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1500
Dune::IndexSet
Index Set Interface base class.
Definition: common/grid.hh:346
Dune::YGridList::end
Iterator end() const
return iterator pointing to the end of the container
Definition: ygrid.hh:927
Dune::Torus::recv
void recv(int rank, void *buffer, int size) const
store a receive request; buffers are received in order; handles also local requests with memcpy
Definition: torus.hh:372
Dune::CommunicationDirection
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
Dune::YaspGrid::GridFamily
YaspGridFamily< dim, Coordinates > GridFamily
the GridFamily of this grid
Definition: yaspgrid.hh:717
Dune::YaspGrid::size
int size(int level, int codim) const
number of entities per level and codim in this process
Definition: yaspgrid.hh:1433
Dune::Overlap_All_Interface
@ Overlap_All_Interface
send overlap, receive all entities
Definition: gridenums.hh:88
Dune::YaspGeometry
The general version that handles all codimensions but 0 and dim.
Definition: yaspgrid.hh:59
Dune::Capabilities::isCartesian::v
static const bool v
Definition: common/capabilities.hh:48
Dune::YGridList::Intersection
type describing an intersection with a neighboring processor
Definition: ygrid.hh:826
Dune::InterfaceType
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
Dune::EquidistantOffsetCoordinates
Container for equidistant coordinates in a YaspGrid with non-trivial origin.
Definition: coordinates.hh:128
Dune::YaspGrid::leafbegin
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1365
Dune::YaspGrid::overlapSize
int overlapSize(int level, int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1407
grid.hh
Different resources needed by all grid implementations.
Dune::YaspGrid::domainSize
const Dune::FieldVector< ctype, dim > & domainSize() const
returns the size of the physical domain
Definition: yaspgrid.hh:1471
Dune
Include standard header files.
Definition: agrid.hh:58
Dune::YaspGrid::Entity
friend class Entity
Definition: yaspgrid.hh:1786
Dune::Capabilities::hasBackupRestoreFacilities
Specialize with 'true' if implementation provides backup and restore facilities. (default=false)
Definition: common/capabilities.hh:121
Dune::VTK::GeometryType
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:180
yaspgridintersection.hh
The YaspIntersection class.
Dune::YaspGrid::LevelIndexSetType
YaspIndexSet< YaspGrid< dim, Coordinates >, false > LevelIndexSetType
Definition: yaspgrid.hh:722
Dune::Capabilities::hasSingleGeometryType::topologyId
static const unsigned int topologyId
Definition: common/capabilities.hh:29
Dune::Entity
Wrapper class for entities.
Definition: common/entity.hh:63
Dune::YaspGrid::GlobalIdSetType
YaspGlobalIdSet< YaspGrid< dim, Coordinates > > GlobalIdSetType
Definition: yaspgrid.hh:724
Dune::YaspGrid::entity
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Definition: yaspgrid.hh:1394
Dune::YaspGrid::ghostSize
int ghostSize(int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1427
Dune::YaspGrid::communicate
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir) const
Definition: yaspgrid.hh:1490