3 #ifndef DUNE_GRID_YASPGRID_HH
4 #define DUNE_GRID_YASPGRID_HH
19 #include <dune/grid/yaspgrid/grids.hh>
21 #include <dune/common/shared_ptr.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/geometry/genericgeometry/topologytypes.hh>
28 #include <dune/geometry/axisalignedcubegeometry.hh>
34 #include <dune/common/parallel/mpicollectivecommunication.hh>
63 template<
int mydim,
int cdim,
class Gr
idImp>
class YaspGeometry;
64 template<
int codim,
int dim,
class Gr
idImp>
class YaspEntity;
65 template<
int codim,
class Gr
idImp>
class YaspEntityPointer;
66 template<
int codim,
class Gr
idImp>
class YaspEntitySeed;
67 template<
int codim, PartitionIteratorType pitype,
class Gr
idImp>
class YaspLevelIterator;
68 template<
class Gr
idImp>
class YaspIntersectionIterator;
69 template<
class Gr
idImp>
class YaspIntersection;
70 template<
class Gr
idImp>
class YaspHierarchicIterator;
71 template<
class Gr
idImp,
bool isLeafIndexSet>
class YaspIndexSet;
72 template<
class Gr
idImp>
class YaspGlobalIdSet;
74 namespace FacadeOptions
77 template<
int dim,
int mydim,
int cdim>
80 static const bool v =
false;
83 template<
int dim,
int mydim,
int cdim>
86 static const bool v =
false;
93 #include <dune/grid/yaspgrid/yaspgridgeometry.hh>
94 #include <dune/grid/yaspgrid/yaspgridentity.hh>
95 #include <dune/grid/yaspgrid/yaspgridintersection.hh>
96 #include <dune/grid/yaspgrid/yaspgridintersectioniterator.hh>
97 #include <dune/grid/yaspgrid/yaspgridhierarchiciterator.hh>
98 #include <dune/grid/yaspgrid/yaspgridentityseed.hh>
99 #include <dune/grid/yaspgrid/yaspgridentitypointer.hh>
100 #include <dune/grid/yaspgrid/yaspgridleveliterator.hh>
101 #include <dune/grid/yaspgrid/yaspgridindexsets.hh>
102 #include <dune/grid/yaspgrid/yaspgrididset.hh>
110 typedef CollectiveCommunication<MPI_Comm>
CCType;
112 typedef CollectiveCommunication<Dune::YaspGrid<dim> >
CCType;
118 YaspGeometry,YaspEntity,
123 YaspIntersectionIterator,
124 YaspIntersectionIterator,
125 YaspHierarchicIterator,
127 YaspIndexSet< const YaspGrid< dim >,
false >,
128 YaspIndexSet< const YaspGrid< dim >,
true >,
129 YaspGlobalIdSet<const YaspGrid<dim> >,
130 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+yaspgrid_codim_bits>,
131 YaspGlobalIdSet<const YaspGrid<dim> >,
132 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+yaspgrid_codim_bits>,
139 template<
int dim,
int codim>
141 template<
class G,
class DataHandle>
144 if (data.contains(dim,codim))
146 DUNE_THROW(
GridError,
"interface communication not implemented");
154 template<
class G,
class DataHandle>
157 if (data.contains(dim,dim))
158 g.template communicateCodim<DataHandle,dim>(data,iftype,dir,level);
165 template<
class G,
class DataHandle>
168 if (data.contains(dim,0))
169 g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
282 DUNE_THROW(
GridError,
"level not existing");
295 static YLoadBalance<dim> lb;
321 for (
int i=0; i<dim; i++) h[i] = L[i]/s[i];
322 for (
int i=0; i<dim; i++) r[i] = 0.5*h[i];
328 for (
int i=0; i<dim; i++)
333 o_overlap[i] = o_interior[i]-overlap;
334 s_overlap[i] = s_interior[i]+2*overlap;
340 int max =
std::min(s[i]-1,o_interior[i]+s_interior[i]-1+overlap);
342 s_overlap[i] = max-min+1;
345 g.
cell_overlap = SubYGrid<dim,ctype>(YGrid<dim,ctype>(o_overlap,s_overlap,h,r));
349 for (
int i=0; i<dim; i++) offset[i] = o_interior[i]-o_overlap[i];
350 g.
cell_interior = SubYGrid<dim,ctype>(o_interior,s_interior,offset,s_overlap,h,r);
357 iTupel o_vertex_global, s_vertex_global;
358 for (
int i=0; i<dim; i++) r[i] = 0.0;
361 for (
int i=0; i<dim; i++) o_vertex_global[i] = g.
cell_global.origin(i);
362 for (
int i=0; i<dim; i++) s_vertex_global[i] = g.
cell_global.size(i)+1;
363 g.
vertex_global = YGrid<dim,ctype>(o_vertex_global,s_vertex_global,h,r);
366 iTupel o_vertex_overlapfront;
367 iTupel s_vertex_overlapfront;
368 for (
int i=0; i<dim; i++) o_vertex_overlapfront[i] = g.
cell_overlap.origin(i);
369 for (
int i=0; i<dim; i++) s_vertex_overlapfront[i] = g.
cell_overlap.size(i)+1;
370 g.
vertex_overlapfront = SubYGrid<dim,ctype>(YGrid<dim,ctype>(o_vertex_overlapfront,s_vertex_overlapfront,h,r));
375 for (
int i=0; i<dim; i++)
383 o_vertex_overlap[i] += 1;
384 s_vertex_overlap[i] -= 1;
390 s_vertex_overlap[i] -= 1;
394 offset[i] = o_vertex_overlap[i]-o_vertex_overlapfront[i];
396 g.
vertex_overlap = SubYGrid<dim,ctype>(o_vertex_overlap,s_vertex_overlap,offset,s_vertex_overlapfront,h,r);
399 iTupel o_vertex_interiorborder;
400 iTupel s_vertex_interiorborder;
401 for (
int i=0; i<dim; i++) o_vertex_interiorborder[i] = g.
cell_interior.origin(i);
402 for (
int i=0; i<dim; i++) s_vertex_interiorborder[i] = g.
cell_interior.size(i)+1;
403 for (
int i=0; i<dim; i++) offset[i] = o_vertex_interiorborder[i]-o_vertex_overlapfront[i];
404 g.
vertex_interiorborder = SubYGrid<dim,ctype>(o_vertex_interiorborder,s_vertex_interiorborder,offset,s_vertex_overlapfront,h,r);
409 for (
int i=0; i<dim; i++)
417 o_vertex_interior[i] += 1;
418 s_vertex_interior[i] -= 1;
424 s_vertex_interior[i] -= 1;
427 offset[i] = o_vertex_interior[i]-o_vertex_overlapfront[i];
429 g.
vertex_interior = SubYGrid<dim,ctype>(o_vertex_interior,s_vertex_interior,offset,s_vertex_overlapfront,h,r);
468 std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
471 std::vector<YGrid<dim,ctype> > send_recvgrid(_torus.neighbors());
472 std::vector<YGrid<dim,ctype> > recv_recvgrid(_torus.neighbors());
473 std::vector<YGrid<dim,ctype> > send_sendgrid(_torus.neighbors());
474 std::vector<YGrid<dim,ctype> > recv_sendgrid(_torus.neighbors());
477 std::vector<mpifriendly_ygrid> mpifriendly_send_recvgrid(_torus.neighbors());
478 std::vector<mpifriendly_ygrid> mpifriendly_recv_recvgrid(_torus.neighbors());
479 std::vector<mpifriendly_ygrid> mpifriendly_send_sendgrid(_torus.neighbors());
480 std::vector<mpifriendly_ygrid> mpifriendly_recv_sendgrid(_torus.neighbors());
484 for (
typename Torus<dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
488 iTupel coord = _torus.coord();
491 for (
int k=0; k<dim; k++) nb[k] += delta[k];
494 for (
int k=0; k<dim; k++)
503 if (nb[k]>=_torus.dims(k))
516 send_sendgrid[i.index()] = sendgrid.move(v);
517 send_recvgrid[i.index()] = recvgrid.move(v);
527 for (
typename Torus<dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
529 mpifriendly_send_sendgrid[i.index()] =
mpifriendly_ygrid(send_sendgrid[i.index()]);
530 _torus.send(i.rank(), &mpifriendly_send_sendgrid[i.index()],
sizeof(
mpifriendly_ygrid));
534 for (
typename Torus<dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
535 _torus.recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()],
sizeof(
mpifriendly_ygrid));
541 for (
typename Torus<dim>::ProcListIterator i=_torus.sendbegin(); i!=_torus.sendend(); ++i)
543 mpifriendly_send_recvgrid[i.index()] =
mpifriendly_ygrid(send_recvgrid[i.index()]);
544 _torus.send(i.rank(), &mpifriendly_send_recvgrid[i.index()],
sizeof(
mpifriendly_ygrid));
548 for (
typename Torus<dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
549 _torus.recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()],
sizeof(
mpifriendly_ygrid));
555 for (
typename Torus<dim>::ProcListIterator i=_torus.recvbegin(); i!=_torus.recvend(); ++i)
560 recv_recvgrid[i.index()] = YGrid<dim,ctype>(yg.
origin,yg.
size,yg.
h,yg.
r);
561 send_intersection.
grid = sendgrid.intersection(recv_recvgrid[i.index()]);
562 send_intersection.
rank = i.rank();
563 send_intersection.
distance = i.distance();
564 if (!send_intersection.
grid.empty()) sendlist.push_front(send_intersection);
567 yg = mpifriendly_recv_sendgrid[i.index()];
568 recv_sendgrid[i.index()] = YGrid<dim,ctype>(yg.
origin,yg.
size,yg.
h,yg.
r);
569 recv_intersection.
grid = recvgrid.intersection(recv_sendgrid[i.index()]);
570 recv_intersection.
rank = i.rank();
571 recv_intersection.
distance = i.distance();
572 if(!recv_intersection.
grid.empty()) recvlist.push_back(recv_intersection);
583 indexsets.push_back( make_shared< YaspIndexSet<
const YaspGrid<dim>,
false > >(*
this,0) );
590 const FieldVector<int, dim> &
size =
begin()->cell_overlap.size();
591 Dune::array<int, dim> sides;
593 for (
int i=0; i<dim; i++)
596 ((
begin()->cell_overlap.origin(i) ==
begin()->cell_global.origin(i))+
597 (
begin()->cell_overlap.origin(i) +
begin()->cell_overlap.size(i)
598 ==
begin()->cell_global.origin(i) +
begin()->cell_global.size(i)));
602 for (
int k=0; k<dim; k++)
605 for (
int l=0; l<dim; l++)
610 nBSegments += sides[k]*offset;
617 typedef bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+yaspgrid_codim_bits>
PersistentIndexType;
630 typedef typename SubYGrid<dim,ctype>::TransformingSubIterator
TSI;
631 typedef typename std::deque<Intersection>::const_iterator
ISIT;
639 _periodic = periodic;
648 array<int,dim> sArray;
649 std::copy(s.begin(), s.end(), sArray.begin());
650 double imbal = _torus.partition(_torus.rank(),o,sArray,o_interior,s_interior);
651 imbal = _torus.global_max(imbal);
658 _levels[0] =
makelevel(0,L,s,periodic,o_interior,s_interior,overlap);
664 Dune::array<int,dim> s,
665 std::bitset<dim> periodic,
670 _periodic = periodic;
674 std::copy(s.begin(), s.end(), this->_s.begin());
680 std::copy(s.begin(), s.end(), s_interior.begin());
682 double imbal = _torus.partition(_torus.rank(),o,s,o_interior,s_interior);
683 imbal = _torus.global_max(imbal);
687 _levels[0] =
makelevel(0,L,_s,periodic,o_interior,s_interior,overlap);
702 Dune::FieldVector<ctype, dim> L,
703 Dune::FieldVector<int, dim> s,
704 Dune::FieldVector<bool, dim> periodic,
int overlap,
706 DUNE_DEPRECATED_MSG("Use the corresponding constructor taking array<
int> and std::bitset")
709 _torus(comm,
tag,s,lb),
713 leafIndexSet_(*
this),
714 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
719 for (
size_t i=0; i<dim; i++)
720 this->_periodic[i] = periodic[i];
741 Dune::FieldVector<int, dim> s,
742 Dune::FieldVector<bool, dim> periodic,
int overlap,
744 DUNE_DEPRECATED_MSG("Use the corresponding constructor taking array<
int> and std::bitset")
746 : ccobj(MPI_COMM_SELF),
747 _torus(MPI_COMM_SELF,
tag,s,lb),
751 leafIndexSet_(*
this),
752 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
757 for (
size_t i=0; i<dim; i++)
758 this->_periodic[i] = periodic[i];
771 Dune::FieldVector<ctype, dim> L,
772 Dune::array<int, dim> s,
773 std::bitset<dim> periodic,
778 _torus(comm,
tag,s,lb),
782 leafIndexSet_(*
this),
783 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
804 Dune::array<int, dim> s,
805 std::bitset<dim> periodic,
809 : ccobj(MPI_COMM_SELF),
810 _torus(MPI_COMM_SELF,
tag,s,lb),
814 leafIndexSet_(*
this),
815 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
832 Dune::array<int, dim> elements)
834 : ccobj(MPI_COMM_SELF),
839 leafIndexSet_(*
this),
843 adaptRefCount(0), adaptActive(
false)
847 std::copy(elements.begin(), elements.end(), _s.begin());
853 std::copy(elements.begin(), elements.end(), s_interior.begin());
855 double imbal = _torus.partition(_torus.rank(),o,elements,o_interior,s_interior);
856 imbal = _torus.global_max(imbal);
860 _levels[0] =
makelevel(0,L,_s,_periodic,o_interior,s_interior,0);
876 return _levels.size()-1;
884 "Coarsening " << -refCount <<
" levels requested!");
887 for (
int k=refCount; k<0; k++)
891 _levels.back() = empty;
896 indexsets.pop_back();
900 for (
int k=0; k<refCount; k++)
907 for (
int i=0; i<dim; i++)
916 for (
int i=0; i<dim; i++)
918 for (
int i=0; i<dim; i++)
922 _levels.push_back(
makelevel(_levels.size(),_LL,s,_periodic,o_interior,s_interior,overlap) );
935 keep_ovlp = keepPhysicalOverlap;
951 assert(adaptActive ==
false);
952 if (e.level() !=
maxLevel())
return false;
953 adaptRefCount =
std::max(adaptRefCount, refCount);
965 return ( e.level() ==
maxLevel() ) ? adaptRefCount : 0;
972 return (adaptRefCount > 0);
979 adaptRefCount =
comm().max(adaptRefCount);
980 return (adaptRefCount < 0);
991 template<
int cd, PartitionIteratorType pitype>
994 return levelbegin<cd,pitype>(level);
998 template<
int cd, PartitionIteratorType pitype>
1001 return levelend<cd,pitype>(level);
1008 return levelbegin<cd,All_Partition>(level);
1015 return levelend<cd,All_Partition>(level);
1019 template<
int cd, PartitionIteratorType pitype>
1022 return levelbegin<cd,pitype>(
maxLevel());
1026 template<
int cd, PartitionIteratorType pitype>
1029 return levelend<cd,pitype>(
maxLevel());
1036 return levelbegin<cd,All_Partition>(
maxLevel());
1043 return levelend<cd,All_Partition>(
maxLevel());
1047 template <
typename Seed>
1048 typename Traits::template Codim<Seed::codimension>::EntityPointer
1051 const int codim = Seed::codimension;
1056 return YaspEntityPointer<codim,GridImp>(
this,g,
1057 TSI(g->cell_overlap, this->getRealImplementation(seed).coord()));
1059 return YaspEntityPointer<codim,GridImp>(
this,g,
1060 TSI(g->vertex_overlap, this->getRealImplementation(seed).coord()));
1062 DUNE_THROW(
GridError,
"YaspEntityPointer: codim not implemented");
1093 int size (
int level,
int codim)
const
1095 return sizes[level][codim];
1107 return (type.isCube()) ? sizes[level][dim-type.dim()] : 0;
1126 template<
class DataHandleImp,
class DataType>
1136 template<
class DataHandleImp,
class DataType>
1146 template<
class DataHandle,
int codim>
1150 if (!data.contains(dim,codim))
return;
1153 typedef typename DataHandle::DataType DataType;
1159 const std::deque<Intersection>* sendlist=0;
1160 const std::deque<Intersection>* recvlist=0;
1167 sendlist = &g->send_cell_interior_overlap;
1168 recvlist = &g->recv_cell_overlap_interior;
1172 sendlist = &g->send_cell_overlap_overlap;
1173 recvlist = &g->recv_cell_overlap_overlap;
1180 sendlist = &g->send_vertex_interiorborder_interiorborder;
1181 recvlist = &g->recv_vertex_interiorborder_interiorborder;
1186 sendlist = &g->send_vertex_interiorborder_overlapfront;
1187 recvlist = &g->recv_vertex_overlapfront_interiorborder;
1191 sendlist = &g->send_vertex_overlap_overlapfront;
1192 recvlist = &g->recv_vertex_overlapfront_overlap;
1196 sendlist = &g->send_vertex_overlapfront_overlapfront;
1197 recvlist = &g->recv_vertex_overlapfront_overlapfront;
1203 std::swap(sendlist,recvlist);
1208 std::vector<int> send_size(sendlist->size(),-1);
1209 std::vector<int> recv_size(recvlist->size(),-1);
1210 std::vector<size_t*> send_sizes(sendlist->size(),
static_cast<size_t*
>(0));
1211 std::vector<size_t*> recv_sizes(recvlist->size(),
static_cast<size_t*
>(0));
1212 if (data.fixedsize(dim,codim))
1216 for (
ISIT is=sendlist->begin(); is!=sendlist->end(); ++is)
1219 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1220 send_size[cnt] = is->grid.totalsize() * data.size(*it);
1224 for (
ISIT is=recvlist->begin(); is!=recvlist->end(); ++is)
1227 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1228 recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1236 for (
ISIT is=sendlist->begin(); is!=sendlist->end(); ++is)
1239 size_t *buf =
new size_t[is->grid.totalsize()];
1240 send_sizes[cnt] = buf;
1243 int i=0;
size_t n=0;
1245 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1247 tsubend(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubend()));
1248 for ( ; it!=tsubend; ++it)
1250 buf[i] = data.size(*it);
1259 torus().send(is->rank,buf,is->grid.totalsize()*
sizeof(size_t));
1265 for (
ISIT is=recvlist->begin(); is!=recvlist->end(); ++is)
1268 size_t *buf =
new size_t[is->grid.totalsize()];
1269 recv_sizes[cnt] = buf;
1272 torus().recv(is->rank,buf,is->grid.totalsize()*
sizeof(size_t));
1281 for (
ISIT is=sendlist->begin(); is!=sendlist->end(); ++is)
1283 delete[] send_sizes[cnt];
1284 send_sizes[cnt] = 0;
1290 for (
ISIT is=recvlist->begin(); is!=recvlist->end(); ++is)
1293 size_t *buf = recv_sizes[cnt];
1297 for (
int i=0; i<is->grid.totalsize(); ++i)
1308 std::vector<DataType*> sends(sendlist->size(),
static_cast<DataType*
>(0));
1310 for (
ISIT is=sendlist->begin(); is!=sendlist->end(); ++is)
1313 DataType *buf =
new DataType[send_size[cnt]];
1319 MessageBuffer<DataType> mb(buf);
1323 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1325 tsubend(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubend()));
1326 for ( ; it!=tsubend; ++it)
1327 data.gather(mb,*it);
1330 torus().send(is->rank,buf,send_size[cnt]*
sizeof(DataType));
1335 std::vector<DataType*> recvs(recvlist->size(),
static_cast<DataType*
>(0));
1337 for (
ISIT is=recvlist->begin(); is!=recvlist->end(); ++is)
1340 DataType *buf =
new DataType[recv_size[cnt]];
1346 torus().recv(is->rank,buf,recv_size[cnt]*
sizeof(DataType));
1355 for (
ISIT is=sendlist->begin(); is!=sendlist->end(); ++is)
1357 delete[] sends[cnt];
1364 for (
ISIT is=recvlist->begin(); is!=recvlist->end(); ++is)
1367 DataType *buf = recvs[cnt];
1370 MessageBuffer<DataType> mb(buf);
1373 if (data.fixedsize(dim,codim))
1376 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1377 size_t n=data.size(*it);
1379 tsubend(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubend()));
1380 for ( ; it!=tsubend; ++it)
1381 data.scatter(mb,*it,n);
1386 size_t *sbuf = recv_sizes[cnt];
1388 it(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubbegin()));
1390 tsubend(YaspLevelIterator<codim,All_Partition,GridImp>(
this,g,is->grid.tsubend()));
1391 for ( ; it!=tsubend; ++it)
1392 data.scatter(mb,*it,sbuf[i++]);
1405 return theglobalidset;
1410 return theglobalidset;
1415 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1416 return *(indexsets[level]);
1421 return leafIndexSet_;
1427 const CollectiveCommunication<MPI_Comm>&
comm ()
const
1434 const CollectiveCommunication<YaspGrid>&
comm ()
const
1446 friend class Dune::YaspIndexSet<const Dune::
YaspGrid<dim>, true >;
1447 friend class Dune::YaspIndexSet<const Dune::YaspGrid<dim>, false >;
1448 friend class Dune::YaspGlobalIdSet<const Dune::YaspGrid<dim> >;
1450 friend class Dune::YaspIntersectionIterator<const Dune::YaspGrid<dim> >;
1451 friend class Dune::YaspIntersection<const Dune::YaspGrid<dim> >;
1452 friend class Dune::YaspEntity<0, dim, const Dune::YaspGrid<dim> >;
1454 template <int codim_, class GridImp_>
1457 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1461 class MessageBuffer {
1464 MessageBuffer (DT *p)
1473 void write (
const Y& data)
1475 dune_static_assert(( is_same<DT,Y>::value ),
"DataType mismatch");
1481 void read (Y& data)
const
1483 dune_static_assert(( is_same<DT,Y>::value ),
"DataType mismatch");
1498 sizes[g->level()][0] = 1;
1499 for (
int i=0; i<dim; ++i)
1500 sizes[g->level()][0] *= g->cell_overlap.size(i);
1505 sizes[g->level()][1] = 0;
1506 for (
int i=0; i<dim; ++i)
1508 int s=g->cell_overlap.size(i)+1;
1509 for (
int j=0; j<dim; ++j)
1511 s *= g->cell_overlap.size(j);
1512 sizes[g->level()][1] += s;
1519 sizes[g->level()][dim-1] = 0;
1520 for (
int i=0; i<dim; ++i)
1522 int s=g->cell_overlap.size(i);
1523 for (
int j=0; j<dim; ++j)
1525 s *= g->cell_overlap.size(j)+1;
1526 sizes[g->level()][dim-1] += s;
1531 sizes[g->level()][dim] = 1;
1532 for (
int i=0; i<dim; ++i)
1533 sizes[g->level()][dim] *= g->vertex_overlapfront.size(i);
1538 template<
int cd, PartitionIteratorType pitype>
1539 YaspLevelIterator<cd,pitype,GridImp> levelbegin (
int level)
const
1541 dune_static_assert( cd == dim || cd == 0 ,
1542 "YaspGrid only supports Entities with codim=dim and codim=0");
1544 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1546 return levelend <cd, pitype> (level);
1550 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->cell_interior.tsubbegin());
1552 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->cell_overlap.tsubbegin());
1557 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_interior.tsubbegin());
1559 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_interiorborder.tsubbegin());
1561 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_overlap.tsubbegin());
1563 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_overlapfront.tsubbegin());
1565 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1569 template<
int cd, PartitionIteratorType pitype>
1570 YaspLevelIterator<cd,pitype,GridImp> levelend (
int level)
const
1572 dune_static_assert( cd == dim || cd == 0 ,
1573 "YaspGrid only supports Entities with codim=dim and codim=0");
1575 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1579 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->cell_interior.tsubend());
1581 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->cell_overlap.tsubend());
1586 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_interior.tsubend());
1588 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_interiorborder.tsubend());
1590 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_overlap.tsubend());
1592 return YaspLevelIterator<cd,pitype,GridImp>(
this,g,g->vertex_overlapfront.tsubend());
1594 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1598 CollectiveCommunication<MPI_Comm> ccobj;
1600 CollectiveCommunication<YaspGrid> ccobj;
1605 std::vector< shared_ptr< YaspIndexSet<const YaspGrid<dim>,
false > > > indexsets;
1606 YaspIndexSet<const YaspGrid<dim>,
true> leafIndexSet_;
1607 YaspGlobalIdSet<const YaspGrid<dim> > theglobalidset;
1611 std::bitset<dim> _periodic;
1612 ReservedVector<YGridLevel,32> _levels;
1614 int sizes[32][dim+1];
1623 inline std::ostream& operator<< (std::ostream& s, YaspGrid<d>& grid)
1625 int rank = grid.torus().rank();
1627 s <<
"[" << rank <<
"]:" <<
" YaspGrid maxlevel=" << grid.maxLevel() << std::endl;
1631 s <<
"[" << rank <<
"]: " << std::endl;
1632 s <<
"[" << rank <<
"]: " <<
"==========================================" << std::endl;
1633 s <<
"[" << rank <<
"]: " <<
"level=" << g->level() << std::endl;
1634 s <<
"[" << rank <<
"]: " <<
"cell_global=" << g->cell_global << std::endl;
1635 s <<
"[" << rank <<
"]: " <<
"cell_overlap=" << g->cell_overlap << std::endl;
1636 s <<
"[" << rank <<
"]: " <<
"cell_interior=" << g->cell_interior << std::endl;
1638 i!=g->send_cell_overlap_overlap.end(); ++i)
1640 s <<
"[" << rank <<
"]: " <<
" s_c_o_o "
1641 << i->rank <<
" " << i->grid << std::endl;
1644 i!=g->recv_cell_overlap_overlap.end(); ++i)
1646 s <<
"[" << rank <<
"]: " <<
" r_c_o_o "
1647 << i->rank <<
" " << i->grid << std::endl;
1650 i!=g->send_cell_interior_overlap.end(); ++i)
1652 s <<
"[" << rank <<
"]: " <<
" s_c_i_o "
1653 << i->rank <<
" " << i->grid << std::endl;
1656 i!=g->recv_cell_overlap_interior.end(); ++i)
1658 s <<
"[" << rank <<
"]: " <<
" r_c_o_i "
1659 << i->rank <<
" " << i->grid << std::endl;
1662 s <<
"[" << rank <<
"]: " <<
"-----------------------------------------------" << std::endl;
1663 s <<
"[" << rank <<
"]: " <<
"vertex_global=" << g->vertex_global << std::endl;
1664 s <<
"[" << rank <<
"]: " <<
"vertex_overlapfront=" << g->vertex_overlapfront << std::endl;
1665 s <<
"[" << rank <<
"]: " <<
"vertex_overlap=" << g->vertex_overlap << std::endl;
1666 s <<
"[" << rank <<
"]: " <<
"vertex_interiorborder=" << g->vertex_interiorborder << std::endl;
1667 s <<
"[" << rank <<
"]: " <<
"vertex_interior=" << g->vertex_interior << std::endl;
1668 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->send_vertex_overlapfront_overlapfront.begin();
1669 i!=g->send_vertex_overlapfront_overlapfront.end(); ++i)
1671 s <<
"[" << rank <<
"]: " <<
" s_v_of_of "
1672 << i->rank <<
" " << i->grid << std::endl;
1674 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->recv_vertex_overlapfront_overlapfront.begin();
1675 i!=g->recv_vertex_overlapfront_overlapfront.end(); ++i)
1677 s <<
"[" << rank <<
"]: " <<
" r_v_of_of "
1678 << i->rank <<
" " << i->grid << std::endl;
1680 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->send_vertex_overlap_overlapfront.begin();
1681 i!=g->send_vertex_overlap_overlapfront.end(); ++i)
1683 s <<
"[" << rank <<
"]: " <<
" s_v_o_of "
1684 << i->rank <<
" " << i->grid << std::endl;
1686 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->recv_vertex_overlapfront_overlap.begin();
1687 i!=g->recv_vertex_overlapfront_overlap.end(); ++i)
1689 s <<
"[" << rank <<
"]: " <<
" r_v_of_o "
1690 << i->rank <<
" " << i->grid << std::endl;
1692 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->send_vertex_interiorborder_interiorborder.begin();
1693 i!=g->send_vertex_interiorborder_interiorborder.end(); ++i)
1695 s <<
"[" << rank <<
"]: " <<
" s_v_ib_ib "
1696 << i->rank <<
" " << i->grid << std::endl;
1698 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->recv_vertex_interiorborder_interiorborder.begin();
1699 i!=g->recv_vertex_interiorborder_interiorborder.end(); ++i)
1701 s <<
"[" << rank <<
"]: " <<
" r_v_ib_ib "
1702 << i->rank <<
" " << i->grid << std::endl;
1704 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->send_vertex_interiorborder_overlapfront.begin();
1705 i!=g->send_vertex_interiorborder_overlapfront.end(); ++i)
1707 s <<
"[" << rank <<
"]: " <<
" s_v_ib_of "
1708 << i->rank <<
" " << i->grid << std::endl;
1710 for (
typename std::deque<
typename YaspGrid<d>::Intersection>::const_iterator i=g->recv_vertex_overlapfront_interiorborder.begin();
1711 i!=g->recv_vertex_overlapfront_interiorborder.end(); ++i)
1713 s <<
"[" << rank <<
"]: " <<
" s_v_of_ib "
1714 << i->rank <<
" " << i->grid << std::endl;
1723 namespace Capabilities
1740 static const bool v =
true;
1741 static const unsigned int topologyId = GenericGeometry :: CubeTopology< dim > :: type :: id ;
1750 static const bool v =
true;
1759 static const bool v =
true;
1768 static const bool v =
true;
1774 static const bool v =
true;
1780 static const bool v =
true;
1789 static const bool v =
true;
1798 static const bool v =
true;
1807 static const bool v =
true;