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