3 #ifndef DUNE_GRID_YASPGRID_HH
4 #define DUNE_GRID_YASPGRID_HH
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>
35 #include <dune/common/parallel/mpicollectivecommunication.hh>
56 template<
int dim,
class Coordinates>
class YaspGrid;
58 template<
int codim,
int dim,
class Gr
idImp>
class YaspEntity;
88 template<
int dim,
class Coordinates>
92 typedef CollectiveCommunication<MPI_Comm>
CCType;
94 typedef CollectiveCommunication<No_Comm>
CCType;
106 YaspIntersectionIterator,
110 YaspIndexSet< const YaspGrid< dim, Coordinates >,
true >,
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>,
122 template<
int dim,
int codim>
123 struct YaspCommunicateMeta {
124 template<
class G,
class DataHandle>
127 if (data.contains(dim,codim))
129 g.template communicateCodim<DataHandle,codim>(data,iftype,dir,level);
131 YaspCommunicateMeta<dim,codim-1>::comm(g,data,iftype,dir,level);
136 struct YaspCommunicateMeta<dim,0> {
137 template<
class G,
class DataHandle>
140 if (data.contains(dim,0))
141 g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
163 template<
int dim,
class Coordinates = Equ
idistantCoordinates<
double, dim> >
165 :
public GridDefaultImplementation<dim,dim,typename Coordinates::ctype,YaspGridFamily<dim, Coordinates> >
168 template<
int, PartitionIteratorType,
typename>
180 typedef typename Coordinates::ctype
ctype;
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;
210 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interior_data;
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;
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;
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;
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;
242 typedef std::array<int, dim> iTupel;
243 typedef FieldVector<ctype, dim> fTupel;
270 return _coarseSize[i] * (1 << l);
277 for (
int i=0; i<dim; ++i)
303 YGridLevelIterator
begin (
int i)
const
306 DUNE_THROW(
GridError,
"level not existing");
311 YGridLevelIterator
end ()
const
331 void makelevel (
const Coordinates& coords, std::bitset<dim> periodic, iTupel o_interior,
int overlap)
333 YGridLevel& g = _levels.back();
338 g.keepOverlap = keep_ovlp;
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();
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();
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();
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();
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();
367 std::array<int,dim> n;
368 std::fill(n.begin(), n.end(), 0);
371 std::bitset<dim> ovlp_low(0ULL);
372 std::bitset<dim> ovlp_up(0ULL);
378 for (
int i=0; i<dim; i++)
382 s_overlap[i] = g.coords.size(i);
387 o_overlap[i] = o_interior[i]-
overlap;
394 if (o_interior[i] - overlap < 0)
398 o_overlap[i] = o_interior[i] -
overlap;
403 if (o_overlap[i] + g.coords.size(i) <
globalSize(i))
408 for (
unsigned int codim = 0; codim < dim + 1; codim++)
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);
425 for (
unsigned int index = 0; index < (1<<dim); index++)
428 std::bitset<dim> r(index);
429 if (r.count() != dim-codim)
433 std::array<int,dim> origin(o_overlap);
434 std::array<int,dim>
size(s_overlap);
438 for (
int i=0; i<dim; i++)
444 for (
int i=0; i<dim; i++)
460 for (
int i=0; i<dim; i++)
482 for (
int i=0; i<dim; i++)
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);
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;
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]);
542 struct mpifriendly_ygrid {
545 std::fill(origin.begin(), origin.end(), 0);
546 std::fill(
size.begin(),
size.end(), 0);
548 mpifriendly_ygrid (
const YGridComponent<Coordinates>& grid)
549 : origin(grid.origin()),
size(grid.
size())
565 std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
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());
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());
587 iTupel coord = _torus.
coord();
588 iTupel delta = i.delta();
590 for (
int k=0; k<dim; k++) nb[k] += delta[k];
592 std::fill(v.begin(), v.end(), 0);
594 for (
int k=0; k<dim; k++)
603 if (nb[k]>=_torus.
dims(k))
616 send_sendgrid[i.index()] = sendgrid.
move(v);
617 send_recvgrid[i.index()] = recvgrid.
move(v);
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));
635 _torus.
recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
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));
649 _torus.
recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
658 Intersection send_intersection;
659 mpifriendly_ygrid yg = mpifriendly_recv_recvgrid[i.index()];
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);
666 Intersection recv_intersection;
667 yg = mpifriendly_recv_sendgrid[i.index()];
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);
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();
692 std::array<int, dim> sides;
694 for (
int i=0; i<dim; i++)
697 ((
begin()->overlap[0].dataBegin()->origin(i) == 0)+
698 (
begin()->overlap[0].dataBegin()->origin(i) +
begin()->overlap[0].dataBegin()->size(i)
703 for (
int k=0; k<dim; k++)
706 for (
int l=0; l<dim; l++)
709 offset *=
begin()->overlap[0].dataBegin()->size(l);
711 nBSegments += sides[k]*offset;
739 std::array<int, dim> s,
740 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
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)
750 "YaspGrid coordinate container template parameter and given constructor values do not match!");
755 std::fill(o.begin(), o.end(), 0);
756 iTupel o_interior(o);
757 iTupel s_interior(s);
762 for (
int i=0; i<dim; i++)
763 if ((s_interior[i] <=
overlap) &&
764 (periodic[i] || (s_interior[i] != s[i])))
768 for (
int i=0; i<dim; i++)
771 iTupel s_overlap(s_interior);
772 for (
int i=0; i<dim; i++)
774 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
776 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
798 Dune::FieldVector<ctype, dim> upperright,
799 std::array<int, dim> s,
800 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
804 : ccobj(
comm), _torus(
comm,tag,s,lb), leafIndexSet_(*this),
805 _L(upperright - lowerleft),
806 _periodic(periodic), _coarseSize(s), _overlap(
overlap),
807 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
811 "YaspGrid coordinate container template parameter and given constructor values do not match!");
816 std::fill(o.begin(), o.end(), 0);
817 iTupel o_interior(o);
818 iTupel s_interior(s);
823 for (
int i=0; i<dim; i++)
824 if ((s_interior[i] <=
overlap) &&
825 (periodic[i] || (s_interior[i] != s[i])))
828 Dune::FieldVector<ctype,dim> extension(upperright);
829 Dune::FieldVector<ctype,dim> h;
830 for (
int i=0; i<dim; i++)
832 extension[i] -= lowerleft[i];
833 h[i] = extension[i] / s[i];
836 iTupel s_overlap(s_interior);
837 for (
int i=0; i<dim; i++)
839 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
841 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
860 YaspGrid (std::array<std::vector<ctype>, dim> coords,
861 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
866 leafIndexSet_(*this), _periodic(periodic), _overlap(
overlap),
867 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
870 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
874 "YaspGrid coordinate container template parameter and given constructor values do not match!");
879 for (
int i=0; i<dim; i++) {
880 _coarseSize[i] = coords[i].size() - 1;
881 _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
885 std::fill(o.begin(), o.end(), 0);
886 iTupel o_interior(o);
887 iTupel s_interior(_coarseSize);
889 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
892 for (
int i=0; i<dim; i++)
893 if ((s_interior[i] <=
overlap) &&
894 (periodic[i] || (s_interior[i] != _coarseSize[i])))
897 std::array<std::vector<ctype>,dim> newcoords;
898 std::array<int, dim> offset(o_interior);
901 for (
int i=0; i<dim; ++i)
904 typename std::vector<ctype>::iterator
begin = coords[i].begin() + o_interior[i];
905 typename std::vector<ctype>::iterator
end = begin + s_interior[i] + 1;
909 if (o_interior[i] -
overlap > 0)
914 if (o_interior[i] + s_interior[i] +
overlap < _coarseSize[i])
918 newcoords[i].resize(end-begin);
919 std::copy(begin, end, newcoords[i].
begin());
923 if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
926 typename std::vector<ctype>::iterator it = coords[i].begin();
928 newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
931 if ((periodic[i]) && (o_interior[i] - overlap <= 0))
936 typename std::vector<ctype>::iterator it = coords[i].end() - 1;
938 newcoords[i].insert(newcoords[i].
begin(), newcoords[i].
front() - *it + *(--it));
960 DUNE_DEPRECATED_MSG(
"This Yaspgrid constructor is deprecated.")
962 Dune::FieldVector<ctype, dim> L,
963 std::array<
int, dim> s,
964 std::bitset<dim> periodic,
967 : ccobj(comm), _torus(comm,tag,s,lb), leafIndexSet_(*this),
968 _L(L), keep_ovlp(true), adaptRefCount(0), adaptActive(false)
970 _periodic = periodic;
976 std::fill(o.begin(), o.end(), 0);
977 iTupel o_interior(o);
978 iTupel s_interior(s);
983 for (
int i=0; i<dim; i++)
986 iTupel s_overlap(s_interior);
987 for (
int i=0; i<dim; i++)
989 if ((o_interior[i] - overlap > 0) || (periodic[i]))
991 if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
997 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1002 makelevel(cc,periodic,o_interior,overlap);
1017 DUNE_DEPRECATED_MSG(
"This Yaspgrid constructor is deprecated.")
1019 std::array<
std::vector<ctype>, dim> coords,
1022 : ccobj(comm), _torus(comm,tag,
Dune::Yasp::
sizeArray<dim>(coords),lb),
1023 leafIndexSet_(*this),
1024 _periodic(
std::bitset<dim>(0)),
1027 adaptRefCount(0), adaptActive(false)
1030 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1031 _periodic = periodic;
1036 for (
int i=0; i<dim; i++) {
1037 _coarseSize[i] = coords[i].size() - 1;
1038 _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
1042 std::fill(o.begin(), o.end(), 0);
1043 iTupel o_interior(o);
1044 iTupel s_interior(_coarseSize);
1046 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
1048 std::array<std::vector<ctype>,dim> newcoords;
1049 std::array<int, dim> offset(o_interior);
1052 for (
int i=0; i<dim; ++i)
1055 typename std::vector<ctype>::iterator
begin = coords[i].begin() + o_interior[i];
1056 typename std::vector<ctype>::iterator
end = begin + s_interior[i] + 1;
1060 if (o_interior[i] - overlap > 0)
1065 if (o_interior[i] + s_interior[i] + overlap < _coarseSize[i])
1069 newcoords[i].resize(end-begin);
1070 std::copy(begin, end, newcoords[i].
begin());
1074 if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
1077 typename std::vector<ctype>::iterator it = coords[i].begin();
1079 newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
1082 if ((periodic[i]) && (o_interior[i] - overlap <= 0))
1087 typename std::vector<ctype>::iterator it = coords[i].end() - 1;
1089 newcoords[i].insert(newcoords[i].
begin(), newcoords[i].
front() - *it + *(--it));
1095 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1100 makelevel(cc,periodic,o_interior,overlap);
1120 YaspGrid (std::array<std::vector<ctype>, dim> coords,
1121 std::bitset<dim> periodic,
1123 CollectiveCommunicationType
comm,
1124 std::array<int,dim> coarseSize,
1126 : ccobj(comm), _torus(comm,tag,coarseSize,lb), leafIndexSet_(*this),
1127 _periodic(periodic), _coarseSize(coarseSize), _overlap(overlap),
1128 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
1132 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1135 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1137 for (
int i=0; i<dim; i++)
1138 _L[i] = coords[i][coords[i].
size() - 1] - coords[i][0];
1142 std::array<int,dim> o;
1143 std::fill(o.begin(), o.end(), 0);
1144 std::array<int,dim> o_interior(o);
1145 std::array<int,dim> s_interior(coarseSize);
1147 _torus.
partition(_torus.
rank(),o,coarseSize,o_interior,s_interior);
1150 std::array<int,dim> offset(o_interior);
1151 for (
int i=0; i<dim; i++)
1152 if ((periodic[i]) || (o_interior[i] > 0))
1158 makelevel(cc,periodic,o_interior,overlap);
1176 return _levels.size()-1;
1184 "Coarsening " << -refCount <<
" levels requested!");
1187 for (
int k=refCount; k<0; k++)
1191 _levels.back() = empty;
1195 indexsets.pop_back();
1199 for (
int k=0; k<refCount; k++)
1202 YGridLevel& cg = _levels[
maxLevel()];
1204 std::bitset<dim> ovlp_low(0ULL), ovlp_up(0ULL);
1205 for (
int i=0; i<dim; i++)
1207 if (cg.overlap[0].dataBegin()->origin(i) > 0 || _periodic[i])
1209 if (cg.overlap[0].dataBegin()->max(i) + 1 <
globalSize(i) || _periodic[i])
1213 Coordinates newcont(cg.coords.refine(ovlp_low, ovlp_up, cg.overlapSize, keep_ovlp));
1215 int overlap = (keep_ovlp) ? 2*cg.overlapSize : cg.overlapSize;
1219 for (
int i=0; i<dim; i++)
1220 o_interior[i] = 2*cg.interior[0].dataBegin()->origin(i);
1223 _levels.resize(_levels.size() + 1);
1224 makelevel(newcont,_periodic,o_interior,overlap);
1236 keep_ovlp = keepPhysicalOverlap;
1252 assert(adaptActive ==
false);
1253 if (e.level() !=
maxLevel())
return false;
1254 adaptRefCount =
std::max(adaptRefCount, refCount);
1266 return ( e.level() ==
maxLevel() ) ? adaptRefCount : 0;
1273 return (adaptRefCount > 0);
1280 adaptRefCount =
comm().max(adaptRefCount);
1281 return (adaptRefCount < 0);
1287 adaptActive =
false;
1292 template<
int cd, PartitionIteratorType pitype>
1295 return levelbegin<cd,pitype>(level);
1299 template<
int cd, PartitionIteratorType pitype>
1302 return levelend<cd,pitype>(level);
1309 return levelbegin<cd,All_Partition>(level);
1316 return levelend<cd,All_Partition>(level);
1320 template<
int cd, PartitionIteratorType pitype>
1323 return levelbegin<cd,pitype>(
maxLevel());
1327 template<
int cd, PartitionIteratorType pitype>
1330 return levelend<cd,pitype>(
maxLevel());
1337 return levelbegin<cd,All_Partition>(
maxLevel());
1344 return levelend<cd,All_Partition>(
maxLevel());
1355 template <
typename Seed>
1356 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.")
1357 typename Traits::template Codim<Seed::codimension>::
EntityPointer
1360 const int codim = Seed::codimension;
1368 template <
typename Seed>
1369 typename Traits::template Codim<Seed::codimension>::Entity
1372 const int codim = Seed::codimension;
1379 return Entity(EntityImp(g,YIterator(g->overlapfront[codim],this->getRealImplementation(seed).coord(),this->
getRealImplementation(seed).offset())));
1385 YGridLevelIterator g =
begin(level);
1386 return g->overlapSize;
1393 return g->overlapSize;
1409 int size (
int level,
int codim)
const
1411 YGridLevelIterator g =
begin(level);
1415 typedef typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator DAI;
1416 for (DAI it = g->overlapfront[codim].dataBegin(); it != g->overlapfront[codim].dataEnd(); ++it)
1417 count += it->totalsize();
1431 return (type.isCube()) ?
size(level,dim-type.dim()) : 0;
1455 template<
class DataHandleImp,
class DataType>
1458 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,level);
1465 template<
class DataHandleImp,
class DataType>
1468 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,this->
maxLevel());
1475 template<
class DataHandle,
int codim>
1479 if (!data.contains(dim,codim))
return;
1482 typedef typename DataHandle::DataType DataType;
1485 YGridLevelIterator g =
begin(level);
1493 sendlist = &g->send_interiorborder_interiorborder[codim];
1494 recvlist = &g->recv_interiorborder_interiorborder[codim];
1498 sendlist = &g->send_interiorborder_overlapfront[codim];
1499 recvlist = &g->recv_overlapfront_interiorborder[codim];
1503 sendlist = &g->send_overlap_overlapfront[codim];
1504 recvlist = &g->recv_overlapfront_overlap[codim];
1508 sendlist = &g->send_overlapfront_overlapfront[codim];
1509 recvlist = &g->recv_overlapfront_overlapfront[codim];
1514 std::swap(sendlist,recvlist);
1519 std::vector<int> send_size(sendlist->
size(),-1);
1520 std::vector<int> recv_size(recvlist->
size(),-1);
1521 std::vector<size_t*> send_sizes(sendlist->
size(),
static_cast<size_t*
>(0));
1522 std::vector<size_t*> recv_sizes(recvlist->
size(),
static_cast<size_t*
>(0));
1527 if (data.fixedsize(dim,codim))
1531 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1535 send_size[cnt] = is->grid.totalsize() * data.size(*it);
1539 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1543 recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1551 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1554 size_t *buf =
new size_t[is->grid.totalsize()];
1555 send_sizes[cnt] = buf;
1558 int i=0;
size_t n=0;
1563 for ( ; it!=itend; ++it)
1565 buf[i] = data.size(*it);
1574 torus().
send(is->rank,buf,is->grid.totalsize()*
sizeof(size_t));
1580 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1583 size_t *buf =
new size_t[is->grid.totalsize()];
1584 recv_sizes[cnt] = buf;
1587 torus().
recv(is->rank,buf,is->grid.totalsize()*
sizeof(size_t));
1596 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1598 delete[] send_sizes[cnt];
1599 send_sizes[cnt] = 0;
1605 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1608 size_t *buf = recv_sizes[cnt];
1612 for (
int i=0; i<is->grid.totalsize(); ++i)
1623 std::vector<DataType*> sends(sendlist->
size(),
static_cast<DataType*
>(0));
1625 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1628 DataType *buf =
new DataType[send_size[cnt]];
1634 MessageBuffer<DataType> mb(buf);
1641 for ( ; it!=itend; ++it)
1642 data.gather(mb,*it);
1645 torus().
send(is->rank,buf,send_size[cnt]*
sizeof(DataType));
1650 std::vector<DataType*> recvs(recvlist->
size(),
static_cast<DataType*
>(0));
1652 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1655 DataType *buf =
new DataType[recv_size[cnt]];
1661 torus().
recv(is->rank,buf,recv_size[cnt]*
sizeof(DataType));
1670 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1672 delete[] sends[cnt];
1679 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1682 DataType *buf = recvs[cnt];
1685 MessageBuffer<DataType> mb(buf);
1688 if (data.fixedsize(dim,codim))
1692 size_t n=data.size(*it);
1695 for ( ; it!=itend; ++it)
1696 data.scatter(mb,*it,n);
1701 size_t *sbuf = recv_sizes[cnt];
1706 for ( ; it!=itend; ++it)
1707 data.scatter(mb,*it,sbuf[i++]);
1720 return theglobalidset;
1725 return theglobalidset;
1730 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1731 return *(indexsets[level]);
1736 return leafIndexSet_;
1741 const CollectiveCommunicationType&
comm ()
const
1761 template <int codim_, class GridImp_>
1764 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1768 class MessageBuffer {
1771 MessageBuffer (DT *p)
1780 void write (
const Y& data)
1782 static_assert(( is_same<DT,Y>::value ),
"DataType mismatch");
1788 void read (Y& data)
const
1790 static_assert(( is_same<DT,Y>::value ),
"DataType mismatch");
1801 template<
int cd, PartitionIteratorType pitype>
1804 YGridLevelIterator g =
begin(level);
1805 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1816 return levelend <cd, pitype> (level);
1818 DUNE_THROW(
GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1822 template<
int cd, PartitionIteratorType pitype>
1823 YaspLevelIterator<cd,pitype,GridImp> levelend (
int level)
const
1825 YGridLevelIterator g =
begin(level);
1826 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1829 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].end());
1831 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].end());
1833 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].end());
1835 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].end());
1837 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1840 CollectiveCommunicationType ccobj;
1842 Torus<CollectiveCommunicationType,dim> _torus;
1844 std::vector< std::shared_ptr< YaspIndexSet<const YaspGrid<dim,Coordinates>,
false > > > indexsets;
1845 YaspIndexSet<const YaspGrid<dim,Coordinates>,
true> leafIndexSet_;
1846 YaspGlobalIdSet<const YaspGrid<dim,Coordinates> > theglobalidset;
1848 Dune::FieldVector<ctype, dim> _L;
1850 std::bitset<dim> _periodic;
1852 ReservedVector<YGridLevel,32> _levels;
1861 template <
int d,
class CC>
1862 std::ostream& operator<< (std::ostream& s, const YaspGrid<d,CC>& grid)
1864 int rank = grid.torus().rank();
1866 s <<
"[" << rank <<
"]:" <<
" YaspGrid maxlevel=" << grid.maxLevel() << std::endl;
1868 s <<
"Printing the torus: " <<std::endl;
1869 s << grid.torus() << std::endl;
1873 s <<
"[" << rank <<
"]: " << std::endl;
1874 s <<
"[" << rank <<
"]: " <<
"==========================================" << std::endl;
1875 s <<
"[" << rank <<
"]: " <<
"level=" << g->level() << std::endl;
1877 for (
int codim = 0; codim < d + 1; ++codim)
1879 s <<
"[" << rank <<
"]: " <<
"overlapfront[" << codim <<
"]: " << g->overlapfront[codim] << std::endl;
1880 s <<
"[" << rank <<
"]: " <<
"overlap[" << codim <<
"]: " << g->overlap[codim] << std::endl;
1881 s <<
"[" << rank <<
"]: " <<
"interiorborder[" << codim <<
"]: " << g->interiorborder[codim] << std::endl;
1882 s <<
"[" << rank <<
"]: " <<
"interior[" << codim <<
"]: " << g->interior[codim] << std::endl;
1885 for (I i=g->send_overlapfront_overlapfront[codim].begin();
1886 i!=g->send_overlapfront_overlapfront[codim].end(); ++i)
1887 s <<
"[" << rank <<
"]: " <<
" s_of_of[" << codim <<
"] to rank "
1888 << i->rank <<
" " << i->grid << std::endl;
1890 for (I i=g->recv_overlapfront_overlapfront[codim].begin();
1891 i!=g->recv_overlapfront_overlapfront[codim].end(); ++i)
1892 s <<
"[" << rank <<
"]: " <<
" r_of_of[" << codim <<
"] to rank "
1893 << i->rank <<
" " << i->grid << std::endl;
1895 for (I i=g->send_overlap_overlapfront[codim].begin();
1896 i!=g->send_overlap_overlapfront[codim].end(); ++i)
1897 s <<
"[" << rank <<
"]: " <<
" s_o_of[" << codim <<
"] to rank "
1898 << i->rank <<
" " << i->grid << std::endl;
1900 for (I i=g->recv_overlapfront_overlap[codim].begin();
1901 i!=g->recv_overlapfront_overlap[codim].end(); ++i)
1902 s <<
"[" << rank <<
"]: " <<
" r_of_o[" << codim <<
"] to rank "
1903 << i->rank <<
" " << i->grid << std::endl;
1905 for (I i=g->send_interiorborder_interiorborder[codim].begin();
1906 i!=g->send_interiorborder_interiorborder[codim].end(); ++i)
1907 s <<
"[" << rank <<
"]: " <<
" s_ib_ib[" << codim <<
"] to rank "
1908 << i->rank <<
" " << i->grid << std::endl;
1910 for (I i=g->recv_interiorborder_interiorborder[codim].begin();
1911 i!=g->recv_interiorborder_interiorborder[codim].end(); ++i)
1912 s <<
"[" << rank <<
"]: " <<
" r_ib_ib[" << codim <<
"] to rank "
1913 << i->rank <<
" " << i->grid << std::endl;
1915 for (I i=g->send_interiorborder_overlapfront[codim].begin();
1916 i!=g->send_interiorborder_overlapfront[codim].end(); ++i)
1917 s <<
"[" << rank <<
"]: " <<
" s_ib_of[" << codim <<
"] to rank "
1918 << i->rank <<
" " << i->grid << std::endl;
1920 for (I i=g->recv_overlapfront_interiorborder[codim].begin();
1921 i!=g->recv_overlapfront_interiorborder[codim].end(); ++i)
1922 s <<
"[" << rank <<
"]: " <<
" r_of_ib[" << codim <<
"] to rank "
1923 << i->rank <<
" " << i->grid << std::endl;
1932 namespace Capabilities
1942 template<
int dim,
class Coordinates>
1945 static const bool v =
true;
1951 template<
int dim,
class Coordinates>
1954 static const bool v =
true;
1955 static const unsigned int topologyId = GenericGeometry :: CubeTopology< dim > :: type :: id ;
1961 template<
int dim,
class Coordinates>
1964 static const bool v =
true;
1970 template<
int dim,
class Coordinates,
int codim>
1973 static const bool v =
true;
1979 template<
int dim,
int codim,
class Coordinates>
1982 static const bool v =
true;
1989 template<
int dim,
class Coordinates>
1990 struct DUNE_DEPRECATED_MSG("Capabilities::
isParallel will be removed after dune-grid-2.4.")
isParallel< YaspGrid<dim, Coordinates> >
1992 static const bool DUNE_DEPRECATED_MSG(
"Capabilities::isParallel will be removed after dune-grid-2.4.") v = true;
1998 template<
int dim, class Coordinates>
2001 static const bool v =
true;
2007 template<
int dim,
class Coordinates>
2010 static const bool v =
true;
type describing an intersection with a neighboring processor
Definition: ygrid.hh:826
void communicateCodim(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1476
CollectiveCommunication< MPI_Comm > CCType
Definition: yaspgrid.hh:92
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:1270
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lbegin(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1307
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
Definition: defaultgridview.hh:223
Wrapper class for pointers to entities.
Definition: common/entitypointer.hh:112
Definition: common/geometry.hh:24
int overlapSize(int level, int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1383
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
interior, border, and overlap entities
Definition: gridenums.hh:137
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1321
const YaspGrid< dim, Coordinates > GridImp
Definition: yaspgrid.hh:678
send interior and border, receive all entities
Definition: gridenums.hh:86
bool getRefineOption() const
Definition: yaspgrid.hh:288
[ provides Dune::Grid ]
Definition: yaspgrid.hh:56
int rank() const
return own rank
Definition: torus.hh:96
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:437
int neighbors() const
return the number of neighbors, which is
Definition: torus.hh:205
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
The YaspIntersection class.
unsigned char uint8_t
Definition: yaspgrid.hh:15
size_t numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: yaspgrid.hh:1441
bool preAdapt()
returns true, if the grid will be coarsened
Definition: yaspgrid.hh:1277
int levelSize(int l, int i) const
return size of the grid (in cells) on level l in direction i
Definition: yaspgrid.hh:268
bool checkIfMonotonous(const Dune::array< std::vector< ctype >, dim > &coords)
Definition: coordinates.hh:361
YaspHierarchicIterator enables iteration over son entities of codim 0.
Definition: yaspgrid.hh:64
Interior interior
PartitionSet for the interior partition.
Definition: partitionset.hh:226
Dune::array< int, d > sizeArray(const Dune::array< std::vector< ct >, d > &v)
Definition: ygrid.hh:26
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:178
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: common/capabilities.hh:47
void exchange() const
exchange messages stored in request buffers; clear request buffers afterwards
Definition: torus.hh:389
Specialization of the PersistentContainer for YaspGrid.
Definition: yaspgrid.hh:67
The general version that handles all codimensions but 0 and dim.
Definition: yaspgrid.hh:57
int max(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:335
Include standard header files.
Definition: agrid.hh:59
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
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:860
The YaspEntityPointer class.
send/receive interior and border entities
Definition: gridenums.hh:85
Definition: defaultgridview.hh:23
static const YLoadBalanceDefault< dim > * defaultLoadbalancer()
Definition: yaspgrid.hh:317
Front front
PartitionSet for the front partition.
Definition: partitionset.hh:235
ReservedVector< YGridLevel, 32 >::const_iterator YGridLevelIterator
Iterator over the grid levels.
Definition: yaspgrid.hh:294
send overlap, receive all entities
Definition: gridenums.hh:88
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: yaspgrid.hh:1429
send all and receive all entities
Definition: gridenums.hh:89
int maxLevel() const
Definition: yaspgrid.hh:1174
level-wise, non-persistent, consecutive indices for YaspGrid
int overlapSize(int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1390
friend class Entity
Definition: yaspgrid.hh:1765
Traits::template Codim< Seed::codimension >::EntityPointer entityPointer(const Seed &seed) const
obtain EntityPointer from EntitySeed.
Definition: yaspgrid.hh:1358
send overlap, receive overlap and front entities
Definition: gridenums.hh:87
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:1250
Iterator over a collection o YGrids A YGrid::Iterator is the heart of an entity in YaspGrid...
Definition: ygrid.hh:590
bool isPeriodic(int i) const
return whether the grid is periodic in direction i
Definition: yaspgrid.hh:283
Container for equidistant coordinates in a YaspGrid with non-trivial origin.
Definition: coordinates.hh:124
ProcListIterator sendend() const
end of send list
Definition: torus.hh:345
facility for writing and reading grids
Definition: common/backuprestore.hh:40
Specialize with 'true' for all codims that a grid implements entities for. (default=false) ...
Definition: common/capabilities.hh:57
only ghost entities
Definition: gridenums.hh:140
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:1466
int ghostSize(int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1403
implements a collection of multiple std::deque Intersections with neighboring processor...
Definition: ygrid.hh:820
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: yaspgrid.hh:1435
CommDataHandleIF describes the features of a data handle for communication in parallel runs using the...
Definition: datahandleif.hh:72
CollectiveCommunication< MPI_Comm > CollectiveCommunicationType
Definition: yaspgrid.hh:182
const int yaspgrid_level_bits
Definition: yaspgrid.hh:50
persistent, globally unique Ids
Definition: yaspgrid.hh:66
Different resources needed by all grid implementations.
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
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:1447
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1335
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
const Traits::LocalIdSet & localIdSet() const
Definition: yaspgrid.hh:1723
const Traits::LeafIndexSet & leafIndexSet() const
Definition: yaspgrid.hh:1734
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1342
all entities
Definition: gridenums.hh:139
Definition: yaspgrid.hh:58
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
returns adaptation mark for given entity
Definition: yaspgrid.hh:1264
void globalRefine(int refCount)
refine the grid refCount times.
Definition: yaspgrid.hh:1180
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:1285
The YaspLevelIterator 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
This file provides the infrastructure for toroidal communication in YaspGrid.
Iterates over entities of one grid level.
Definition: yaspgrid.hh:61
int ghostSize(int level, int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1397
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:1423
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
YaspIndexSet< YaspGrid< dim, Coordinates >, true > LeafIndexSetType
Definition: yaspgrid.hh:727
A pointer to a YaspGrid::Entity.
Definition: yaspgrid.hh:59
reverse communication direction
Definition: gridenums.hh:170
Id Set Interface.
Definition: common/grid.hh:362
static const bool v
Definition: common/capabilities.hh:50
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
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Definition: yaspgrid.hh:1370
Overlap overlap
PartitionSet for the overlap partition.
Definition: partitionset.hh:232
static const bool v
Definition: common/capabilities.hh:118
Specialize with 'true' if implementation guarantees conforming level grids. (default=false) ...
Definition: common/capabilities.hh:98
Intersection of a mesh entities of codimension 0 ("elements") with a "neighboring" element or with th...
Definition: albertagrid/dgfparser.hh:26
const iTupel & dims() const
return dimensions of torus
Definition: torus.hh:114
the YaspEntity class and its specializations
Index Set Interface base class.
Definition: common/grid.hh:361
const int yaspgrid_dim_bits
Definition: yaspgrid.hh:49
void refineOptions(bool keepPhysicalOverlap)
set options for refinement
Definition: yaspgrid.hh:1234
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:16
Implement the default load balance strategy of yaspgrid.
Definition: partitioning.hh:34
YaspIndexSet< YaspGrid< dim, Coordinates >, false > LevelIndexSetType
Definition: yaspgrid.hh:726
The YaspGeometry class and its specializations.
const Traits::GlobalIdSet & globalIdSet() const
Definition: yaspgrid.hh:1718
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
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false) ...
Definition: common/capabilities.hh:107
const Traits::LevelIndexSet & levelIndexSet(int level) const
Definition: yaspgrid.hh:1728
Container for equidistant coordinates in a YaspGrid.
Definition: coordinates.hh:26
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lend(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1314
The YaspIntersectionIterator class.
ProcListIterator recvend() const
last process in receive list
Definition: torus.hh:357
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
static const bool v
Definition: common/capabilities.hh:28
interior and border entities
Definition: gridenums.hh:136
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1328
static const unsigned int topologyId
Definition: common/capabilities.hh:31
int size(int level, int codim) const
number of entities per level and codim in this process
Definition: yaspgrid.hh:1409
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
static const bool v
Definition: common/capabilities.hh:59
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:1300
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
This provides a YGrid, the elemental component of the yaspgrid implementation.
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1456
Types for GridView.
Definition: common/grid.hh:420
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lbegin(int level) const
one past the end on this level
Definition: yaspgrid.hh:1293
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
Specialize with 'true' if implementation supports parallelism. (default=false)
Definition: common/capabilities.hh:68
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
only interior entities
Definition: gridenums.hh:135
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:797
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:1741
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 grid model.
Definition: common/grid.hh:1343
iTupel globalSize() const
return number of cells on finest level on all processors
Definition: yaspgrid.hh:262
The YaspEntitySeed 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.
Wrapper class for entities.
Definition: common/entity.hh:61
YaspGlobalIdSet< YaspGrid< dim, Coordinates > > GlobalIdSetType
Definition: yaspgrid.hh:728
implements a collection of YGridComponents which form a codimension Entities of given codimension c n...
Definition: ygrid.hh:547