3 #ifndef DUNE_GRID_YASPGRID_HH
4 #define DUNE_GRID_YASPGRID_HH
10 #include <type_traits>
22 #include <dune/common/hybridutilities.hh>
23 #include <dune/common/power.hh>
24 #include <dune/common/bigunsignedint.hh>
25 #include <dune/common/typetraits.hh>
26 #include <dune/common/reservedvector.hh>
27 #include <dune/common/parallel/communication.hh>
28 #include <dune/common/parallel/mpihelper.hh>
29 #include <dune/common/deprecated.hh>
30 #include <dune/geometry/axisalignedcubegeometry.hh>
31 #include <dune/geometry/type.hh>
37 #include <dune/common/parallel/mpicommunication.hh>
58 template<
int dim,
class Coordinates>
class YaspGrid;
59 template<
int mydim,
int cdim,
class Gr
idImp>
class YaspGeometry;
60 template<
int codim,
int dim,
class Gr
idImp>
class YaspEntity;
62 template<
int codim, PartitionIteratorType pitype,
class Gr
idImp>
class YaspLevelIterator;
66 template<
class Gr
idImp,
bool isLeafIndexSet>
class YaspIndexSet;
88 template<
int dim,
class Coordinates>
92 typedef CollectiveCommunication<MPI_Comm>
CCType;
94 typedef CollectiveCommunication<No_Comm>
CCType;
111 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
113 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
121 template<
int dim,
int codim>
122 struct YaspCommunicateMeta {
123 template<
class G,
class DataHandle>
126 if (data.contains(dim,codim))
128 g.template communicateCodim<DataHandle,codim>(data,iftype,dir,level);
130 YaspCommunicateMeta<dim,codim-1>::comm(g,data,iftype,dir,level);
135 struct YaspCommunicateMeta<dim,0> {
136 template<
class G,
class DataHandle>
139 if (data.contains(dim,0))
140 g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
162 template<
int dim,
class Coordinates = Equ
idistantCoordinates<
double, dim> >
167 template<
int, PartitionIteratorType,
typename>
179 typedef typename Coordinates::ctype
ctype;
202 std::array<YGrid, dim+1> overlapfront;
203 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlapfront_data;
204 std::array<YGrid, dim+1>
overlap;
205 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlap_data;
206 std::array<YGrid, dim+1> interiorborder;
207 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interiorborder_data;
209 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interior_data;
211 std::array<YGridList<Coordinates>,dim+1> send_overlapfront_overlapfront;
212 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlapfront_overlapfront_data;
213 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlapfront;
214 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlapfront_data;
216 std::array<YGridList<Coordinates>,dim+1> send_overlap_overlapfront;
217 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlap_overlapfront_data;
218 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlap;
219 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlap_data;
221 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_interiorborder;
222 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_interiorborder_data;
223 std::array<YGridList<Coordinates>,dim+1> recv_interiorborder_interiorborder;
224 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_interiorborder_interiorborder_data;
226 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_overlapfront;
227 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_overlapfront_data;
228 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_interiorborder;
229 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_interiorborder_data;
241 typedef std::array<int, dim> iTupel;
242 typedef FieldVector<ctype, dim> fTupel;
269 return _coarseSize[i] * (1 << l);
276 for (
int i=0; i<dim; ++i)
305 DUNE_THROW(
GridError,
"level not existing");
330 void makelevel (
const Coordinates& coords, std::bitset<dim> periodic, iTupel o_interior,
int overlap)
332 YGridLevel& g = _levels.back();
337 g.keepOverlap = keep_ovlp;
340 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlapfront_it = g.overlapfront_data.begin();
341 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlap_it = g.overlap_data.begin();
342 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interiorborder_it = g.interiorborder_data.begin();
343 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interior_it = g.interior_data.begin();
345 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
346 send_overlapfront_overlapfront_it = g.send_overlapfront_overlapfront_data.begin();
347 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
348 recv_overlapfront_overlapfront_it = g.recv_overlapfront_overlapfront_data.begin();
350 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
351 send_overlap_overlapfront_it = g.send_overlap_overlapfront_data.begin();
352 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
353 recv_overlapfront_overlap_it = g.recv_overlapfront_overlap_data.begin();
355 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
356 send_interiorborder_interiorborder_it = g.send_interiorborder_interiorborder_data.begin();
357 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
358 recv_interiorborder_interiorborder_it = g.recv_interiorborder_interiorborder_data.begin();
360 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
361 send_interiorborder_overlapfront_it = g.send_interiorborder_overlapfront_data.begin();
362 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
363 recv_overlapfront_interiorborder_it = g.recv_overlapfront_interiorborder_data.begin();
366 std::array<int,dim> n;
367 std::fill(n.begin(), n.end(), 0);
370 std::bitset<dim> ovlp_low(0ULL);
371 std::bitset<dim> ovlp_up(0ULL);
377 for (
int i=0; i<dim; i++)
381 s_overlap[i] = g.coords.size(i);
386 o_overlap[i] = o_interior[i]-
overlap;
393 if (o_interior[i] -
overlap < 0)
397 o_overlap[i] = o_interior[i] -
overlap;
402 if (o_overlap[i] + g.coords.size(i) <
globalSize(i))
407 for (
unsigned int codim = 0; codim < dim + 1; codim++)
410 g.overlapfront[codim].setBegin(overlapfront_it);
411 g.overlap[codim].setBegin(overlap_it);
412 g.interiorborder[codim].setBegin(interiorborder_it);
413 g.interior[codim].setBegin(interior_it);
414 g.send_overlapfront_overlapfront[codim].setBegin(send_overlapfront_overlapfront_it);
415 g.recv_overlapfront_overlapfront[codim].setBegin(recv_overlapfront_overlapfront_it);
416 g.send_overlap_overlapfront[codim].setBegin(send_overlap_overlapfront_it);
417 g.recv_overlapfront_overlap[codim].setBegin(recv_overlapfront_overlap_it);
418 g.send_interiorborder_interiorborder[codim].setBegin(send_interiorborder_interiorborder_it);
419 g.recv_interiorborder_interiorborder[codim].setBegin(recv_interiorborder_interiorborder_it);
420 g.send_interiorborder_overlapfront[codim].setBegin(send_interiorborder_overlapfront_it);
421 g.recv_overlapfront_interiorborder[codim].setBegin(recv_overlapfront_interiorborder_it);
424 for (
unsigned int index = 0; index < (1<<dim); index++)
427 std::bitset<dim> r(index);
428 if (r.count() != dim-codim)
432 std::array<int,dim> origin(o_overlap);
433 std::array<int,dim>
size(s_overlap);
437 for (
int i=0; i<dim; i++)
443 for (
int i=0; i<dim; i++)
459 for (
int i=0; i<dim; i++)
481 for (
int i=0; i<dim; i++)
496 intersections(*overlapfront_it,*overlapfront_it,*send_overlapfront_overlapfront_it, *recv_overlapfront_overlapfront_it);
497 intersections(*overlap_it,*overlapfront_it,*send_overlap_overlapfront_it, *recv_overlapfront_overlap_it);
498 intersections(*interiorborder_it,*interiorborder_it,*send_interiorborder_interiorborder_it,*recv_interiorborder_interiorborder_it);
499 intersections(*interiorborder_it,*overlapfront_it,*send_interiorborder_overlapfront_it,*recv_overlapfront_interiorborder_it);
506 ++send_overlapfront_overlapfront_it;
507 ++recv_overlapfront_overlapfront_it;
508 ++send_overlap_overlapfront_it;
509 ++recv_overlapfront_overlap_it;
510 ++send_interiorborder_interiorborder_it;
511 ++recv_interiorborder_interiorborder_it;
512 ++send_interiorborder_overlapfront_it;
513 ++recv_overlapfront_interiorborder_it;
517 g.overlapfront[codim].finalize(overlapfront_it);
518 g.overlap[codim].finalize(overlap_it);
519 g.interiorborder[codim].finalize(interiorborder_it);
520 g.interior[codim].finalize(interior_it);
521 g.send_overlapfront_overlapfront[codim].finalize(send_overlapfront_overlapfront_it,g.overlapfront[codim]);
522 g.recv_overlapfront_overlapfront[codim].finalize(recv_overlapfront_overlapfront_it,g.overlapfront[codim]);
523 g.send_overlap_overlapfront[codim].finalize(send_overlap_overlapfront_it,g.overlapfront[codim]);
524 g.recv_overlapfront_overlap[codim].finalize(recv_overlapfront_overlap_it,g.overlapfront[codim]);
525 g.send_interiorborder_interiorborder[codim].finalize(send_interiorborder_interiorborder_it,g.overlapfront[codim]);
526 g.recv_interiorborder_interiorborder[codim].finalize(recv_interiorborder_interiorborder_it,g.overlapfront[codim]);
527 g.send_interiorborder_overlapfront[codim].finalize(send_interiorborder_overlapfront_it,g.overlapfront[codim]);
528 g.recv_overlapfront_interiorborder[codim].finalize(recv_overlapfront_interiorborder_it,g.overlapfront[codim]);
541 struct mpifriendly_ygrid {
544 std::fill(origin.begin(), origin.end(), 0);
545 std::fill(
size.begin(),
size.end(), 0);
547 mpifriendly_ygrid (
const YGridComponent<Coordinates>& grid)
548 : origin(grid.origin()),
size(grid.
size())
564 std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
569 std::vector<YGridComponent<Coordinates> > send_recvgrid(_torus.
neighbors());
570 std::vector<YGridComponent<Coordinates> > recv_recvgrid(_torus.
neighbors());
571 std::vector<YGridComponent<Coordinates> > send_sendgrid(_torus.
neighbors());
572 std::vector<YGridComponent<Coordinates> > recv_sendgrid(_torus.
neighbors());
575 std::vector<mpifriendly_ygrid> mpifriendly_send_recvgrid(_torus.
neighbors());
576 std::vector<mpifriendly_ygrid> mpifriendly_recv_recvgrid(_torus.
neighbors());
577 std::vector<mpifriendly_ygrid> mpifriendly_send_sendgrid(_torus.
neighbors());
578 std::vector<mpifriendly_ygrid> mpifriendly_recv_sendgrid(_torus.
neighbors());
586 iTupel coord = _torus.
coord();
587 iTupel delta = i.delta();
589 for (
int k=0; k<dim; k++) nb[k] += delta[k];
591 std::fill(v.begin(), v.end(), 0);
593 for (
int k=0; k<dim; k++)
602 if (nb[k]>=_torus.
dims(k))
615 send_sendgrid[i.index()] = sendgrid.
move(v);
616 send_recvgrid[i.index()] = recvgrid.
move(v);
628 mpifriendly_send_sendgrid[i.index()] = mpifriendly_ygrid(send_sendgrid[i.index()]);
629 _torus.
send(i.rank(), &mpifriendly_send_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
634 _torus.
recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
642 mpifriendly_send_recvgrid[i.index()] = mpifriendly_ygrid(send_recvgrid[i.index()]);
643 _torus.
send(i.rank(), &mpifriendly_send_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
648 _torus.
recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
658 mpifriendly_ygrid yg = mpifriendly_recv_recvgrid[i.index()];
660 send_intersection.grid = sendgrid.
intersection(recv_recvgrid[i.index()]);
661 send_intersection.rank = i.rank();
662 send_intersection.distance = i.distance();
663 if (!send_intersection.grid.empty()) sendlist.push_front(send_intersection);
666 yg = mpifriendly_recv_sendgrid[i.index()];
668 recv_intersection.grid = recvgrid.
intersection(recv_sendgrid[i.index()]);
669 recv_intersection.rank = i.rank();
670 recv_intersection.distance = i.distance();
671 if(!recv_intersection.grid.empty()) recvlist.push_back(recv_intersection);
688 std::array<int, dim> sides;
690 for (
int i=0; i<dim; i++)
693 ((
begin()->overlap[0].dataBegin()->origin(i) == 0)+
694 (
begin()->overlap[0].dataBegin()->origin(i) +
begin()->overlap[0].dataBegin()->size(i)
699 for (
int k=0; k<dim; k++)
702 for (
int l=0; l<dim; l++)
705 offset *=
begin()->overlap[0].dataBegin()->size(l);
707 nBSegments += sides[k]*offset;
734 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
739 , leafIndexSet_(*this)
740 , _periodic(periodic)
749 for (std::size_t i=0; i<dim; i++)
750 _coarseSize[i] = coordinates.size(i);
753 _torus = decltype(_torus)(
comm,tag,_coarseSize,lb);
756 std::fill(o.begin(), o.end(), 0);
757 iTupel o_interior(o);
758 iTupel s_interior(_coarseSize);
760 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
766 for (std::size_t i=0; i<dim; i++)
767 _L[i] = coordinates.size(i) * coordinates.meshsize(i,0);
772 for (
int i=0; i<dim; i++)
773 _L[i] = coordinates.coordinate(i,_coarseSize[i]) - coordinates.coordinate(i,0);
778 int mysteryFactor = (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value) ? 1 : 2;
781 for (
int i=0; i<dim; i++)
784 int toosmall = (s_interior[i] <= mysteryFactor *
overlap) &&
785 (periodic[i] || (s_interior[i] != _coarseSize[i]));
788 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
790 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
791 " Note that this also holds for DOFs on subdomain boundaries."
792 " Increase grid elements or decrease overlap accordingly.");
799 iTupel s_overlap(s_interior);
800 for (
int i=0; i<dim; i++)
802 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
804 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
808 FieldVector<ctype,dim> upperRightWithOverlap;
809 for (
int i=0; i<dim; i++)
810 upperRightWithOverlap[i] = coordinates.coordinate(i,0) + coordinates.meshsize(i,0) * s_overlap[i];
823 Dune::FieldVector<ctype,dim> lowerleft;
824 for (
int i=0; i<dim; i++)
825 lowerleft[i] =
id(coordinates).
origin(i);
837 std::array<std::vector<ctype>,dim> newCoords;
838 std::array<int, dim> offset(o_interior);
841 for (
int i=0; i<dim; ++i)
844 std::size_t
begin = o_interior[i];
845 std::size_t
end =
begin + s_interior[i] + 1;
849 if (o_interior[i] -
overlap > 0)
854 if (o_interior[i] + s_interior[i] +
overlap < _coarseSize[i])
859 auto newCoordsIt = newCoords[i].begin();
860 for (std::size_t j=
begin; j<
end; j++)
862 *newCoordsIt = coordinates.coordinate(i, j);
868 if ((periodic[i]) && (o_interior[i] + s_interior[i] +
overlap >= _coarseSize[i]))
872 newCoords[i].push_back(newCoords[i].back() - coordinates.coordinate(i,j) + coordinates.coordinate(i,j+1));
875 if ((periodic[i]) && (o_interior[i] -
overlap <= 0))
880 std::size_t reverseCounter = coordinates.size(i);
881 for (
int j=0; j<
overlap; ++j, --reverseCounter)
882 newCoords[i].insert(newCoords[i].
begin(), newCoords[i].
front()
883 - coordinates.coordinate(i,reverseCounter) + coordinates.coordinate(i,reverseCounter-1));
905 std::array<int, dim> s,
906 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
910 : ccobj(
comm), _torus(
comm,tag,s,lb), leafIndexSet_(*this),
911 _L(L), _periodic(periodic), _coarseSize(s), _overlap(
overlap),
912 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
916 "YaspGrid coordinate container template parameter and given constructor values do not match!");
921 std::fill(o.begin(), o.end(), 0);
922 iTupel o_interior(o);
923 iTupel s_interior(s);
929 for (
int i=0; i<dim; i++)
932 int toosmall = (s_interior[i] / 2 <=
overlap) &&
933 (periodic[i] || (s_interior[i] != s[i]));
936 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
938 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
939 " Note that this also holds for DOFs on subdomain boundaries."
940 " Increase grid elements or decrease overlap accordingly.");
944 iTupel s_overlap(s_interior);
945 for (
int i=0; i<dim; i++)
947 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
949 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
953 FieldVector<ctype,dim> upperRightWithOverlap;
955 for (
int i=0; i<dim; i++)
956 upperRightWithOverlap[i] = (L[i] / s[i]) * s_overlap[i];
977 Dune::FieldVector<ctype, dim> upperright,
978 std::array<int, dim> s,
979 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
983 : ccobj(
comm), _torus(
comm,tag,s,lb), leafIndexSet_(*this),
984 _L(upperright - lowerleft),
985 _periodic(periodic), _coarseSize(s), _overlap(
overlap),
986 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
990 "YaspGrid coordinate container template parameter and given constructor values do not match!");
995 std::fill(o.begin(), o.end(), 0);
996 iTupel o_interior(o);
997 iTupel s_interior(s);
1003 for (
int i=0; i<dim; i++)
1006 int toosmall = (s_interior[i] / 2 <=
overlap) &&
1007 (periodic[i] || (s_interior[i] != s[i]));
1010 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
1012 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1013 " Note that this also holds for DOFs on subdomain boundaries."
1014 " Increase grid elements or decrease overlap accordingly.");
1018 iTupel s_overlap(s_interior);
1019 for (
int i=0; i<dim; i++)
1021 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
1023 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
1027 FieldVector<ctype,dim> upperRightWithOverlap;
1028 for (
int i=0; i<dim; i++)
1029 upperRightWithOverlap[i] = lowerleft[i]
1030 + s_overlap[i] * (upperright[i]-lowerleft[i]) / s[i];
1048 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
1053 leafIndexSet_(*this), _periodic(periodic), _overlap(
overlap),
1054 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
1057 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1061 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1066 for (
int i=0; i<dim; i++) {
1067 _coarseSize[i] = coords[i].size() - 1;
1068 _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
1072 std::fill(o.begin(), o.end(), 0);
1073 iTupel o_interior(o);
1074 iTupel s_interior(_coarseSize);
1076 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
1080 for (
int i=0; i<dim; i++)
1083 int toosmall = (s_interior[i] / 2 <=
overlap) &&
1084 (periodic[i] || (s_interior[i] != _coarseSize[i]));
1087 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
1089 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1090 " Note that this also holds for DOFs on subdomain boundaries."
1091 " Increase grid elements or decrease overlap accordingly.");
1096 std::array<std::vector<ctype>,dim> newcoords;
1097 std::array<int, dim> offset(o_interior);
1100 for (
int i=0; i<dim; ++i)
1103 typename std::vector<ctype>::iterator
begin = coords[i].begin() + o_interior[i];
1104 typename std::vector<ctype>::iterator
end =
begin + s_interior[i] + 1;
1108 if (o_interior[i] -
overlap > 0)
1113 if (o_interior[i] + s_interior[i] +
overlap < _coarseSize[i])
1122 if ((periodic[i]) && (o_interior[i] + s_interior[i] +
overlap >= _coarseSize[i]))
1125 typename std::vector<ctype>::iterator it = coords[i].begin();
1127 newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
1130 if ((periodic[i]) && (o_interior[i] -
overlap <= 0))
1135 typename std::vector<ctype>::iterator it = coords[i].end() - 1;
1137 newcoords[i].insert(newcoords[i].
begin(), newcoords[i].
front() - *it + *(--it));
1165 std::bitset<dim> periodic,
1168 std::array<int,dim> coarseSize,
1170 : ccobj(
comm), _torus(
comm,tag,coarseSize,lb), leafIndexSet_(*this),
1171 _periodic(periodic), _coarseSize(coarseSize), _overlap(
overlap),
1172 keep_ovlp(true), adaptRefCount(0), adaptActive(false)
1176 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1179 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1181 for (
int i=0; i<dim; i++)
1182 _L[i] = coords[i][coords[i].
size() - 1] - coords[i][0];
1186 std::array<int,dim> o;
1187 std::fill(o.begin(), o.end(), 0);
1188 std::array<int,dim> o_interior(o);
1189 std::array<int,dim> s_interior(coarseSize);
1191 _torus.
partition(_torus.
rank(),o,coarseSize,o_interior,s_interior);
1194 std::array<int,dim> offset(o_interior);
1195 for (
int i=0; i<dim; i++)
1196 if ((periodic[i]) || (o_interior[i] > 0))
1220 return _levels.size()-1;
1228 "Coarsening " << -refCount <<
" levels requested!");
1231 for (
int k=refCount; k<0; k++)
1235 _levels.back() = empty;
1239 indexsets.pop_back();
1243 for (
int k=0; k<refCount; k++)
1246 YGridLevel& cg = _levels[
maxLevel()];
1248 std::bitset<dim> ovlp_low(0ULL), ovlp_up(0ULL);
1249 for (
int i=0; i<dim; i++)
1251 if (cg.overlap[0].dataBegin()->origin(i) > 0 || _periodic[i])
1253 if (cg.overlap[0].dataBegin()->max(i) + 1 <
globalSize(i) || _periodic[i])
1257 Coordinates newcont(cg.coords.refine(ovlp_low, ovlp_up, cg.overlapSize, keep_ovlp));
1259 int overlap = (keep_ovlp) ? 2*cg.overlapSize : cg.overlapSize;
1263 for (
int i=0; i<dim; i++)
1264 o_interior[i] = 2*cg.interior[0].dataBegin()->origin(i);
1267 _levels.resize(_levels.size() + 1);
1280 keep_ovlp = keepPhysicalOverlap;
1294 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e )
1296 assert(adaptActive ==
false);
1297 if (e.level() !=
maxLevel())
return false;
1298 adaptRefCount =
std::max(adaptRefCount, refCount);
1308 int getMark (
const typename Traits::template Codim<0>::Entity &e )
const
1310 return ( e.level() ==
maxLevel() ) ? adaptRefCount : 0;
1317 return (adaptRefCount > 0);
1324 adaptRefCount =
comm().max(adaptRefCount);
1325 return (adaptRefCount < 0);
1331 adaptActive =
false;
1336 template<
int cd, PartitionIteratorType pitype>
1337 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator
lbegin (
int level)
const
1339 return levelbegin<cd,pitype>(level);
1343 template<
int cd, PartitionIteratorType pitype>
1344 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator
lend (
int level)
const
1346 return levelend<cd,pitype>(level);
1351 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator
lbegin (
int level)
const
1353 return levelbegin<cd,All_Partition>(level);
1358 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator
lend (
int level)
const
1360 return levelend<cd,All_Partition>(level);
1364 template<
int cd, PartitionIteratorType pitype>
1365 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator
leafbegin ()
const
1367 return levelbegin<cd,pitype>(
maxLevel());
1371 template<
int cd, PartitionIteratorType pitype>
1372 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator
leafend ()
const
1374 return levelend<cd,pitype>(
maxLevel());
1379 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator
leafbegin ()
const
1381 return levelbegin<cd,All_Partition>(
maxLevel());
1386 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator
leafend ()
const
1388 return levelend<cd,All_Partition>(
maxLevel());
1392 template <
typename Seed>
1393 typename Traits::template Codim<Seed::codimension>::Entity
1396 const int codim = Seed::codimension;
1399 typedef typename Traits::template Codim<Seed::codimension>::Entity
Entity;
1403 return Entity(EntityImp(g,YIterator(g->overlapfront[codim],seed.impl().coord(),seed.impl().offset())));
1410 return g->overlapSize;
1417 return g->overlapSize;
1433 int size (
int level,
int codim)
const
1439 typedef typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator DAI;
1440 for (DAI it = g->overlapfront[codim].dataBegin(); it != g->overlapfront[codim].dataEnd(); ++it)
1441 count += it->totalsize();
1455 return (type.isCube()) ?
size(level,dim-type.dim()) : 0;
1479 template<
class DataHandleImp,
class DataType>
1482 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,level);
1489 template<
class DataHandleImp,
class DataType>
1492 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,this->
maxLevel());
1499 template<
class DataHandle,
int codim>
1503 if (!data.contains(dim,codim))
return;
1506 typedef typename DataHandle::DataType DataType;
1517 sendlist = &g->send_interiorborder_interiorborder[codim];
1518 recvlist = &g->recv_interiorborder_interiorborder[codim];
1522 sendlist = &g->send_interiorborder_overlapfront[codim];
1523 recvlist = &g->recv_overlapfront_interiorborder[codim];
1527 sendlist = &g->send_overlap_overlapfront[codim];
1528 recvlist = &g->recv_overlapfront_overlap[codim];
1532 sendlist = &g->send_overlapfront_overlapfront[codim];
1533 recvlist = &g->recv_overlapfront_overlapfront[codim];
1538 std::swap(sendlist,recvlist);
1543 std::vector<int> send_size(sendlist->
size(),-1);
1544 std::vector<int> recv_size(recvlist->
size(),-1);
1545 std::vector<size_t*> send_sizes(sendlist->
size(),
static_cast<size_t*
>(0));
1546 std::vector<size_t*> recv_sizes(recvlist->
size(),
static_cast<size_t*
>(0));
1551 if (data.fixedSize(dim,codim))
1555 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1557 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1559 send_size[cnt] = is->grid.totalsize() * data.size(*it);
1563 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1565 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1567 recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1575 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1578 size_t *buf =
new size_t[is->grid.totalsize()];
1579 send_sizes[cnt] = buf;
1582 int i=0;
size_t n=0;
1583 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1585 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1587 for ( ; it!=itend; ++it)
1589 buf[i] = data.size(*it);
1598 torus().
send(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1604 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1607 size_t *buf =
new size_t[is->grid.totalsize()];
1608 recv_sizes[cnt] = buf;
1611 torus().
recv(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1620 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1622 delete[] send_sizes[cnt];
1623 send_sizes[cnt] = 0;
1629 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1632 size_t *buf = recv_sizes[cnt];
1636 for (
int i=0; i<is->grid.totalsize(); ++i)
1647 std::vector<DataType*> sends(sendlist->
size(),
static_cast<DataType*
>(0));
1649 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1652 DataType *buf =
new DataType[send_size[cnt]];
1658 MessageBuffer<DataType> mb(buf);
1661 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1663 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1665 for ( ; it!=itend; ++it)
1666 data.gather(mb,*it);
1669 torus().
send(is->rank,buf,send_size[cnt]*
sizeof(DataType));
1674 std::vector<DataType*> recvs(recvlist->
size(),
static_cast<DataType*
>(0));
1676 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1679 DataType *buf =
new DataType[recv_size[cnt]];
1685 torus().
recv(is->rank,buf,recv_size[cnt]*
sizeof(DataType));
1694 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1696 delete[] sends[cnt];
1703 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1706 DataType *buf = recvs[cnt];
1709 MessageBuffer<DataType> mb(buf);
1712 if (data.fixedSize(dim,codim))
1714 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1716 size_t n=data.size(*it);
1717 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1719 for ( ; it!=itend; ++it)
1720 data.scatter(mb,*it,n);
1725 size_t *sbuf = recv_sizes[cnt];
1726 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1728 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1730 for ( ; it!=itend; ++it)
1731 data.scatter(mb,*it,sbuf[i++]);
1744 return theglobalidset;
1749 return theglobalidset;
1754 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1755 return *(indexsets[level]);
1760 return leafIndexSet_;
1783 friend class
Dune::YaspEntity<0, dim, const Dune::YaspGrid<dim, Coordinates> >;
1785 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1789 class MessageBuffer {
1792 MessageBuffer (DT *p)
1801 void write (
const Y& data)
1803 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1809 void read (Y& data)
const
1811 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1822 template<
int cd, PartitionIteratorType pitype>
1826 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1837 return levelend <cd, pitype> (level);
1839 DUNE_THROW(
GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1843 template<
int cd, PartitionIteratorType pitype>
1844 YaspLevelIterator<cd,pitype,GridImp> levelend (
int level)
const
1847 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1850 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].end());
1852 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].end());
1854 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].end());
1856 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].end());
1858 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1863 Torus<CollectiveCommunicationType,dim> _torus;
1865 std::vector< std::shared_ptr< YaspIndexSet<const YaspGrid<dim,Coordinates>,
false > > > indexsets;
1866 YaspIndexSet<const YaspGrid<dim,Coordinates>,
true> leafIndexSet_;
1867 YaspGlobalIdSet<const YaspGrid<dim,Coordinates> > theglobalidset;
1869 Dune::FieldVector<ctype, dim> _L;
1871 std::bitset<dim> _periodic;
1873 ReservedVector<YGridLevel,32> _levels;
1882 template <
int d,
class CC>
1887 s <<
"[" << rank <<
"]:" <<
" YaspGrid maxlevel=" << grid.
maxLevel() << std::endl;
1889 s <<
"Printing the torus: " <<std::endl;
1890 s << grid.
torus() << std::endl;
1894 s <<
"[" << rank <<
"]: " << std::endl;
1895 s <<
"[" << rank <<
"]: " <<
"==========================================" << std::endl;
1896 s <<
"[" << rank <<
"]: " <<
"level=" << g->level() << std::endl;
1898 for (
int codim = 0; codim < d + 1; ++codim)
1900 s <<
"[" << rank <<
"]: " <<
"overlapfront[" << codim <<
"]: " << g->overlapfront[codim] << std::endl;
1901 s <<
"[" << rank <<
"]: " <<
"overlap[" << codim <<
"]: " << g->overlap[codim] << std::endl;
1902 s <<
"[" << rank <<
"]: " <<
"interiorborder[" << codim <<
"]: " << g->interiorborder[codim] << std::endl;
1903 s <<
"[" << rank <<
"]: " <<
"interior[" << codim <<
"]: " << g->interior[codim] << std::endl;
1906 for (I i=g->send_overlapfront_overlapfront[codim].begin();
1907 i!=g->send_overlapfront_overlapfront[codim].end(); ++i)
1908 s <<
"[" << rank <<
"]: " <<
" s_of_of[" << codim <<
"] to rank "
1909 << i->rank <<
" " << i->grid << std::endl;
1911 for (I i=g->recv_overlapfront_overlapfront[codim].begin();
1912 i!=g->recv_overlapfront_overlapfront[codim].end(); ++i)
1913 s <<
"[" << rank <<
"]: " <<
" r_of_of[" << codim <<
"] to rank "
1914 << i->rank <<
" " << i->grid << std::endl;
1916 for (I i=g->send_overlap_overlapfront[codim].begin();
1917 i!=g->send_overlap_overlapfront[codim].end(); ++i)
1918 s <<
"[" << rank <<
"]: " <<
" s_o_of[" << codim <<
"] to rank "
1919 << i->rank <<
" " << i->grid << std::endl;
1921 for (I i=g->recv_overlapfront_overlap[codim].begin();
1922 i!=g->recv_overlapfront_overlap[codim].end(); ++i)
1923 s <<
"[" << rank <<
"]: " <<
" r_of_o[" << codim <<
"] to rank "
1924 << i->rank <<
" " << i->grid << std::endl;
1926 for (I i=g->send_interiorborder_interiorborder[codim].begin();
1927 i!=g->send_interiorborder_interiorborder[codim].end(); ++i)
1928 s <<
"[" << rank <<
"]: " <<
" s_ib_ib[" << codim <<
"] to rank "
1929 << i->rank <<
" " << i->grid << std::endl;
1931 for (I i=g->recv_interiorborder_interiorborder[codim].begin();
1932 i!=g->recv_interiorborder_interiorborder[codim].end(); ++i)
1933 s <<
"[" << rank <<
"]: " <<
" r_ib_ib[" << codim <<
"] to rank "
1934 << i->rank <<
" " << i->grid << std::endl;
1936 for (I i=g->send_interiorborder_overlapfront[codim].begin();
1937 i!=g->send_interiorborder_overlapfront[codim].end(); ++i)
1938 s <<
"[" << rank <<
"]: " <<
" s_ib_of[" << codim <<
"] to rank "
1939 << i->rank <<
" " << i->grid << std::endl;
1941 for (I i=g->recv_overlapfront_interiorborder[codim].begin();
1942 i!=g->recv_overlapfront_interiorborder[codim].end(); ++i)
1943 s <<
"[" << rank <<
"]: " <<
" r_of_ib[" << codim <<
"] to rank "
1944 << i->rank <<
" " << i->grid << std::endl;
1953 namespace Capabilities
1963 template<
int dim,
class Coordinates>
1966 static const bool v =
true;
1972 template<
int dim,
class Coordinates>
1975 static const bool v =
true;
1976 static const unsigned int topologyId = GeometryTypes::cube(dim).id();
1982 template<
int dim,
class Coordinates>
1985 static const bool v =
true;
1991 template<
int dim,
class Coordinates,
int codim>
1994 static const bool v =
true;
2001 template<
int dim,
class Coordinates,
int codim>
2004 static const bool v =
true;
2010 template<
int dim,
int codim,
class Coordinates>
2013 static const bool v =
true;
2019 template<
int dim,
class Coordinates>
2022 static const bool v =
true;
2028 template<
int dim,
class Coordinates>
2031 static const bool v =
true;
This provides container classes for the coordinates to be used in YaspGrid Upon implementation of the...
Specialization of the StructuredGridFactory class for YaspGrid.
This file provides the infrastructure for toroidal communication in YaspGrid.
the YaspEntity class and its specializations
The YaspEntitySeed class.
The YaspGeometry class and its specializations.
level-wise, non-persistent, consecutive indices for YaspGrid
The YaspIntersection class.
The YaspIntersectionIterator class.
The YaspLevelIterator class.
Specialization of the PersistentContainer for YaspGrid.
This provides a YGrid, the elemental component of the yaspgrid implementation.
Describes the parallel communication interface class for MessageBuffers and DataHandles.
Provides base classes for index and id sets.
unsigned char uint8_t
Definition: yaspgrid.hh:16
std::ostream & operator<<(std::ostream &out, const PartitionType &type)
write a PartitionType to a stream
Definition: gridenums.hh:70
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:168
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:84
@ All_Partition
all entities
Definition: gridenums.hh:139
@ Interior_Partition
only interior entities
Definition: gridenums.hh:135
@ InteriorBorder_Partition
interior and border entities
Definition: gridenums.hh:136
@ Overlap_Partition
interior, border, and overlap entities
Definition: gridenums.hh:137
@ Ghost_Partition
only ghost entities
Definition: gridenums.hh:140
@ BackwardCommunication
reverse communication direction
Definition: gridenums.hh:170
@ InteriorBorder_All_Interface
send interior and border, receive all entities
Definition: gridenums.hh:86
@ All_All_Interface
send all and receive all entities
Definition: gridenums.hh:89
@ Overlap_All_Interface
send overlap, receive all entities
Definition: gridenums.hh:88
@ Overlap_OverlapFront_Interface
send overlap, receive overlap and front entities
Definition: gridenums.hh:87
@ InteriorBorder_InteriorBorder_Interface
send/receive interior and border entities
Definition: gridenums.hh:85
Include standard header files.
Definition: agrid.hh:59
const int yaspgrid_level_bits
Definition: yaspgrid.hh:52
const int yaspgrid_dim_bits
Definition: yaspgrid.hh:51
int max(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:335
constexpr Overlap overlap
PartitionSet for the overlap partition.
Definition: partitionset.hh:276
constexpr Front front
PartitionSet for the front partition.
Definition: partitionset.hh:279
constexpr Interior interior
PartitionSet for the interior partition.
Definition: partitionset.hh:270
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:180
bool checkIfMonotonous(const std::array< std::vector< ctype >, dim > &coords)
Definition: coordinates.hh:370
std::array< int, d > sizeArray(const std::array< std::vector< ct >, d > &v)
Definition: ygrid.hh:27
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:162
facility for writing and reading grids
Definition: common/backuprestore.hh:41
Specialize with 'true' for if the codimension 0 entity of the grid has only one possible geometry typ...
Definition: common/capabilities.hh:25
static const bool v
Definition: common/capabilities.hh:26
static const unsigned int topologyId
Definition: common/capabilities.hh:29
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: common/capabilities.hh:46
static const bool v
Definition: common/capabilities.hh:48
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:56
static const bool v
Definition: common/capabilities.hh:57
specialize with 'true' for all codims that a grid provides an iterator for (default=false)
Definition: common/capabilities.hh:72
static const bool v
Definition: common/capabilities.hh:73
specialize with 'true' for all codims that a grid can communicate data on (default=false)
Definition: common/capabilities.hh:95
static const bool v
Definition: common/capabilities.hh:96
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:104
static const bool v
Definition: common/capabilities.hh:105
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:113
static const bool v
Definition: common/capabilities.hh:114
Specialize with 'true' if implementation provides backup and restore facilities. (default=false)
Definition: common/capabilities.hh:122
static const bool v
Definition: common/capabilities.hh:123
CommDataHandleIF describes the features of a data handle for communication in parallel runs using the...
Definition: datahandleif.hh:76
Definition: defaultgridview.hh:24
Definition: defaultgridview.hh:206
Wrapper class for entities.
Definition: common/entity.hh:64
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:18
Definition: common/grid.hh:855
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:1027
Index Set Interface base class.
Definition: indexidset.hh:76
Id Set Interface.
Definition: indexidset.hh:441
GridFamily::Traits::CollectiveCommunication CollectiveCommunication
A type that is a model of Dune::CollectiveCommunication. It provides a portable way for collective co...
Definition: common/grid.hh:519
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:1065
[ provides Dune::Grid ]
Definition: yaspgrid.hh:165
YaspGridFamily< dim, Coordinates > GridFamily
the GridFamily of this grid
Definition: yaspgrid.hh:717
YaspIndexSet< YaspGrid< dim, Coordinates >, false > LevelIndexSetType
Definition: yaspgrid.hh:722
int ghostSize(int level, int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1421
YaspGrid(std::array< std::vector< ctype >, dim > coords, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Standard constructor for a tensorproduct YaspGrid.
Definition: yaspgrid.hh:1047
void init()
Definition: yaspgrid.hh:679
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: yaspgrid.hh:1453
const Traits::LocalIdSet & localIdSet() const
Definition: yaspgrid.hh:1747
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1386
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Definition: yaspgrid.hh:1394
bool getRefineOption() const
Definition: yaspgrid.hh:287
void globalRefine(int refCount)
refine the grid refCount times.
Definition: yaspgrid.hh:1224
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
returns adaptation mark for given entity
Definition: yaspgrid.hh:1308
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lbegin(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1351
void boundarysegmentssize()
Definition: yaspgrid.hh:685
YaspIndexSet< YaspGrid< dim, Coordinates >, true > LeafIndexSetType
Definition: yaspgrid.hh:723
int ghostSize(int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1427
int overlapSize(int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1414
int globalSize(int i) const
return number of cells on finest level in given direction on all processors
Definition: yaspgrid.hh:255
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1379
const Dune::FieldVector< ctype, dim > & domainSize() const
returns the size of the physical domain
Definition: yaspgrid.hh:1471
void postAdapt()
clean up some markers
Definition: yaspgrid.hh:1329
YaspGridFamily< dim, Coordinates >::Traits Traits
Definition: yaspgrid.hh:719
YGridLevelIterator end() const
return iterator pointing to one past the finest level
Definition: yaspgrid.hh:310
const CollectiveCommunicationType & comm() const
return a collective communication object
Definition: yaspgrid.hh:1765
friend class Entity
Definition: yaspgrid.hh:1786
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: yaspgrid.hh:1459
YaspGrid(Dune::FieldVector< ctype, dim > lowerleft, Dune::FieldVector< ctype, dim > upperright, std::array< int, dim > s, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:976
int maxLevel() const
Definition: yaspgrid.hh:1218
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1372
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1480
void intersections(const YGridComponent< Coordinates > &sendgrid, const YGridComponent< Coordinates > &recvgrid, std::deque< Intersection > &sendlist, std::deque< Intersection > &recvlist)
Construct list of intersections with neighboring processors.
Definition: yaspgrid.hh:563
static const YLoadBalanceDefault< dim > * defaultLoadbalancer()
Definition: yaspgrid.hh:316
bool preAdapt()
returns true, if the grid will be coarsened
Definition: yaspgrid.hh:1321
iTupel levelSize(int l) const
return size vector of the grid (in cells) on level l
Definition: yaspgrid.hh:273
bool mark(int refCount, const typename Traits::template Codim< 0 >::Entity &e)
Marks an entity to be refined/coarsened in a subsequent adapt.
Definition: yaspgrid.hh:1294
int overlapSize(int level, int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1407
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir) const
Definition: yaspgrid.hh:1490
int size(int codim) const
number of leaf entities per codim in this process
Definition: yaspgrid.hh:1447
const Traits::GlobalIdSet & globalIdSet() const
Definition: yaspgrid.hh:1742
bool isPeriodic(int i) const
return whether the grid is periodic in direction i
Definition: yaspgrid.hh:282
void refineOptions(bool keepPhysicalOverlap)
set options for refinement
Definition: yaspgrid.hh:1278
YaspGrid(Dune::FieldVector< ctype, dim > L, std::array< int, dim > s, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:904
int size(int level, int codim) const
number of entities per level and codim in this process
Definition: yaspgrid.hh:1433
YaspGlobalIdSet< YaspGrid< dim, Coordinates > > GlobalIdSetType
Definition: yaspgrid.hh:724
const Torus< CollectiveCommunicationType, dim > & torus() const
return reference to torus
Definition: yaspgrid.hh:249
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lend(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1358
bool adapt()
map adapt to global refine
Definition: yaspgrid.hh:1314
ReservedVector< YGridLevel, 32 >::const_iterator YGridLevelIterator
Iterator over the grid levels.
Definition: yaspgrid.hh:293
size_t numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: yaspgrid.hh:1465
void communicateCodim(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1500
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lend(int level) const
Iterator to one past the last entity of given codim on level for partition type.
Definition: yaspgrid.hh:1344
CollectiveCommunication< MPI_Comm > CollectiveCommunicationType
Definition: yaspgrid.hh:181
void makelevel(const Coordinates &coords, std::bitset< dim > periodic, iTupel o_interior, int overlap)
Make a new YGridLevel structure.
Definition: yaspgrid.hh:330
bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim > PersistentIndexType
Definition: yaspgrid.hh:714
const Traits::LeafIndexSet & leafIndexSet() const
Definition: yaspgrid.hh:1758
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1365
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lbegin(int level) const
one past the end on this level
Definition: yaspgrid.hh:1337
iTupel globalSize() const
return number of cells on finest level on all processors
Definition: yaspgrid.hh:261
const YaspGrid< dim, Coordinates > GridImp
Definition: yaspgrid.hh:677
YGridLevelIterator begin(int i) const
return iterator pointing to given level
Definition: yaspgrid.hh:302
Coordinates::ctype ctype
Type used for coordinates.
Definition: yaspgrid.hh:179
YaspGrid(const Coordinates &coordinates, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:733
const Traits::LevelIndexSet & levelIndexSet(int level) const
Definition: yaspgrid.hh:1752
int levelSize(int l, int i) const
return size of the grid (in cells) on level l in direction i
Definition: yaspgrid.hh:267
YGridLevelIterator begin() const
return iterator pointing to coarsest level
Definition: yaspgrid.hh:296
The general version that handles all codimensions but 0 and dim.
Definition: yaspgridgeometry.hh:29
Definition: yaspgridentity.hh:262
Describes the minimal information necessary to create a fully functional YaspEntity.
Definition: yaspgridentityseed.hh:16
Iterates over entities of one grid level.
Definition: yaspgridleveliterator.hh:17
YaspIntersectionIterator enables iteration over intersections with neighboring codim 0 entities.
Definition: yaspgridintersectioniterator.hh:20
YaspIntersection provides data about intersection with neighboring codim 0 entities.
Definition: yaspgridintersection.hh:20
YaspHierarchicIterator enables iteration over son entities of codim 0.
Definition: yaspgridhierarchiciterator.hh:18
Implementation of Level- and LeafIndexSets for YaspGrid.
Definition: yaspgridindexsets.hh:23
persistent, globally unique Ids
Definition: yaspgrididset.hh:23
Definition: yaspgridpersistentcontainer.hh:33
Definition: yaspgrid.hh:90
GridTraits< dim, dim, Dune::YaspGrid< dim, Coordinates >, YaspGeometry, YaspEntity, YaspLevelIterator, YaspIntersection, YaspIntersection, YaspIntersectionIterator, YaspIntersectionIterator, YaspHierarchicIterator, YaspLevelIterator, YaspIndexSet< const YaspGrid< dim, Coordinates >, false >, YaspIndexSet< const YaspGrid< dim, Coordinates >, true >, YaspGlobalIdSet< const YaspGrid< dim, Coordinates > >, bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim >, YaspGlobalIdSet< const YaspGrid< dim, Coordinates > >, bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim >, CCType, DefaultLevelGridViewTraits, DefaultLeafGridViewTraits, YaspEntitySeed > Traits
Definition: yaspgrid.hh:117
CollectiveCommunication< MPI_Comm > CCType
Definition: yaspgrid.hh:92
Container for equidistant coordinates in a YaspGrid.
Definition: coordinates.hh:27
Container for equidistant coordinates in a YaspGrid with non-trivial origin.
Definition: coordinates.hh:129
ct origin(int d) const
Definition: coordinates.hh:183
Coordinate container for a tensor product YaspGrid.
Definition: coordinates.hh:243
a base class for the yaspgrid partitioning strategy The name might be irritating. It will probably ch...
Definition: partitioning.hh:24
Implement the default load balance strategy of yaspgrid.
Definition: partitioning.hh:35
int neighbors() const
return the number of neighbors, which is
Definition: torus.hh:201
ProcListIterator sendbegin() const
first process in send list
Definition: torus.hh:335
int rank() const
return own rank
Definition: torus.hh:92
ProcListIterator recvbegin() const
first process in receive list
Definition: torus.hh:347
void recv(int rank, void *buffer, int size) const
store a receive request; buffers are received in order; handles also local requests with memcpy
Definition: torus.hh:372
const iTupel & dims() const
return dimensions of torus
Definition: torus.hh:110
ProcListIterator recvend() const
last process in receive list
Definition: torus.hh:353
void send(int rank, void *buffer, int size) const
store a send request; buffers are sent in order; handles also local requests with memcpy
Definition: torus.hh:359
double partition(int rank, iTupel origin_in, iTupel size_in, iTupel &origin_out, iTupel &size_out) const
partition the given grid onto the torus and return the piece of the process with given rank; returns ...
Definition: torus.hh:237
iTupel coord() const
return own coordinates
Definition: torus.hh:98
void exchange() const
exchange messages stored in request buffers; clear request buffers afterwards
Definition: torus.hh:385
ProcListIterator sendend() const
end of send list
Definition: torus.hh:341
YGridComponent< Coordinates > move(iTupel v) const
return grid moved by the vector v
Definition: ygrid.hh:261
YGridComponent< Coordinates > intersection(const YGridComponent< Coordinates > &r) const
Return YGridComponent of supergrid of self which is the intersection of self and another YGridCompone...
Definition: ygrid.hh:269
implements a collection of YGridComponents which form a codimension Entities of given codimension c n...
Definition: ygrid.hh:549
Iterator over a collection o YGrids A YGrid::Iterator is the heart of an entity in YaspGrid.
Definition: ygrid.hh:592
implements a collection of multiple std::deque<Intersection> Intersections with neighboring processor...
Definition: ygrid.hh:821
int size() const
return the size of the container, this is the sum of the sizes of all deques
Definition: ygrid.hh:951
Iterator end() const
return iterator pointing to the end of the container
Definition: ygrid.hh:927
Iterator begin() const
return iterator pointing to the begin of the container
Definition: ygrid.hh:921
type describing an intersection with a neighboring processor
Definition: ygrid.hh:827
A set of traits classes to store static information about grid implementation.
Different resources needed by all grid implementations.