3#ifndef DUNE_GRID_YASPGRID_HH
4#define DUNE_GRID_YASPGRID_HH
16#include <dune/common/hybridutilities.hh>
17#include <dune/common/power.hh>
18#include <dune/common/bigunsignedint.hh>
19#include <dune/common/typetraits.hh>
20#include <dune/common/reservedvector.hh>
21#include <dune/common/parallel/communication.hh>
22#include <dune/common/parallel/mpihelper.hh>
23#include <dune/geometry/axisalignedcubegeometry.hh>
24#include <dune/geometry/type.hh>
30#include <dune/common/parallel/mpicommunication.hh>
51 template<
int dim,
class Coordinates>
class YaspGrid;
52 template<
int mydim,
int cdim,
class Gr
idImp>
class YaspGeometry;
53 template<
int codim,
int dim,
class Gr
idImp>
class YaspEntity;
55 template<
int codim, PartitionIteratorType pitype,
class Gr
idImp>
class YaspLevelIterator;
59 template<
class Gr
idImp,
bool isLeafIndexSet>
class YaspIndexSet;
87 template<
int dim,
class Coordinates>
106 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
108 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
116 template<
int dim,
int codim>
117 struct YaspCommunicateMeta {
118 template<
class G,
class DataHandle>
121 if (data.contains(dim,codim))
123 g.template communicateCodim<DataHandle,codim>(data,iftype,dir,level);
125 YaspCommunicateMeta<dim,codim-1>::comm(g,data,iftype,dir,level);
130 struct YaspCommunicateMeta<dim,0> {
131 template<
class G,
class DataHandle>
132 static void comm (
const G& g, DataHandle& data, InterfaceType iftype, CommunicationDirection dir,
int level)
134 if (data.contains(dim,0))
135 g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
157 template<
int dim,
class Coordinates = Equ
idistantCoordinates<
double, dim> >
162 template<
int, PartitionIteratorType,
typename>
172 typedef typename Coordinates::ctype
ctype;
191 std::array<YGrid, dim+1> overlapfront;
192 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlapfront_data;
193 std::array<YGrid, dim+1> overlap;
194 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> overlap_data;
195 std::array<YGrid, dim+1> interiorborder;
196 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interiorborder_data;
197 std::array<YGrid, dim+1> interior;
198 std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power> interior_data;
200 std::array<YGridList<Coordinates>,dim+1> send_overlapfront_overlapfront;
201 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlapfront_overlapfront_data;
202 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlapfront;
203 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlapfront_data;
205 std::array<YGridList<Coordinates>,dim+1> send_overlap_overlapfront;
206 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_overlap_overlapfront_data;
207 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlap;
208 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_overlap_data;
210 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_interiorborder;
211 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_interiorborder_data;
212 std::array<YGridList<Coordinates>,dim+1> recv_interiorborder_interiorborder;
213 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_interiorborder_interiorborder_data;
215 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_overlapfront;
216 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> send_interiorborder_overlapfront_data;
217 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_interiorborder;
218 std::array<std::deque<Intersection>, StaticPower<2,dim>::power> recv_overlapfront_interiorborder_data;
230 typedef std::array<int, dim> iTupel;
231 typedef FieldVector<ctype, dim> fTupel;
258 return _coarseSize[i] * (1 << l);
265 for (
int i=0; i<dim; ++i)
294 DUNE_THROW(
GridError,
"level not existing");
319 void makelevel (
const Coordinates& coords, std::bitset<dim> periodic, iTupel o_interior,
int overlap)
321 YGridLevel& g = _levels.back();
322 g.overlapSize = overlap;
326 g.keepOverlap = keep_ovlp;
329 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlapfront_it = g.overlapfront_data.begin();
330 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator overlap_it = g.overlap_data.begin();
331 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interiorborder_it = g.interiorborder_data.begin();
332 typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator interior_it = g.interior_data.begin();
334 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
335 send_overlapfront_overlapfront_it = g.send_overlapfront_overlapfront_data.begin();
336 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
337 recv_overlapfront_overlapfront_it = g.recv_overlapfront_overlapfront_data.begin();
339 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
340 send_overlap_overlapfront_it = g.send_overlap_overlapfront_data.begin();
341 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
342 recv_overlapfront_overlap_it = g.recv_overlapfront_overlap_data.begin();
344 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
345 send_interiorborder_interiorborder_it = g.send_interiorborder_interiorborder_data.begin();
346 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
347 recv_interiorborder_interiorborder_it = g.recv_interiorborder_interiorborder_data.begin();
349 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
350 send_interiorborder_overlapfront_it = g.send_interiorborder_overlapfront_data.begin();
351 typename std::array<std::deque<Intersection>, StaticPower<2,dim>::power>::iterator
352 recv_overlapfront_interiorborder_it = g.recv_overlapfront_interiorborder_data.begin();
355 std::array<int,dim> n;
356 std::fill(n.begin(), n.end(), 0);
359 std::bitset<dim> ovlp_low(0ULL);
360 std::bitset<dim> ovlp_up(0ULL);
366 for (
int i=0; i<dim; i++)
370 s_overlap[i] = g.coords.size(i);
375 o_overlap[i] = o_interior[i]-overlap;
382 if (o_interior[i] - overlap < 0)
386 o_overlap[i] = o_interior[i] - overlap;
391 if (o_overlap[i] + g.coords.size(i) <
globalSize(i))
396 for (
unsigned int codim = 0; codim < dim + 1; codim++)
399 g.overlapfront[codim].setBegin(overlapfront_it);
400 g.overlap[codim].setBegin(overlap_it);
401 g.interiorborder[codim].setBegin(interiorborder_it);
402 g.interior[codim].setBegin(interior_it);
403 g.send_overlapfront_overlapfront[codim].setBegin(send_overlapfront_overlapfront_it);
404 g.recv_overlapfront_overlapfront[codim].setBegin(recv_overlapfront_overlapfront_it);
405 g.send_overlap_overlapfront[codim].setBegin(send_overlap_overlapfront_it);
406 g.recv_overlapfront_overlap[codim].setBegin(recv_overlapfront_overlap_it);
407 g.send_interiorborder_interiorborder[codim].setBegin(send_interiorborder_interiorborder_it);
408 g.recv_interiorborder_interiorborder[codim].setBegin(recv_interiorborder_interiorborder_it);
409 g.send_interiorborder_overlapfront[codim].setBegin(send_interiorborder_overlapfront_it);
410 g.recv_overlapfront_interiorborder[codim].setBegin(recv_overlapfront_interiorborder_it);
413 for (
unsigned int index = 0; index < (1<<dim); index++)
416 std::bitset<dim> r(index);
417 if (r.count() != dim-codim)
421 std::array<int,dim> origin(o_overlap);
422 std::array<int,dim>
size(s_overlap);
426 for (
int i=0; i<dim; i++)
432 for (
int i=0; i<dim; i++)
448 for (
int i=0; i<dim; i++)
452 origin[i] += overlap;
470 for (
int i=0; i<dim; i++)
485 intersections(*overlapfront_it,*overlapfront_it,*send_overlapfront_overlapfront_it, *recv_overlapfront_overlapfront_it);
486 intersections(*overlap_it,*overlapfront_it,*send_overlap_overlapfront_it, *recv_overlapfront_overlap_it);
487 intersections(*interiorborder_it,*interiorborder_it,*send_interiorborder_interiorborder_it,*recv_interiorborder_interiorborder_it);
488 intersections(*interiorborder_it,*overlapfront_it,*send_interiorborder_overlapfront_it,*recv_overlapfront_interiorborder_it);
495 ++send_overlapfront_overlapfront_it;
496 ++recv_overlapfront_overlapfront_it;
497 ++send_overlap_overlapfront_it;
498 ++recv_overlapfront_overlap_it;
499 ++send_interiorborder_interiorborder_it;
500 ++recv_interiorborder_interiorborder_it;
501 ++send_interiorborder_overlapfront_it;
502 ++recv_overlapfront_interiorborder_it;
506 g.overlapfront[codim].finalize(overlapfront_it);
507 g.overlap[codim].finalize(overlap_it);
508 g.interiorborder[codim].finalize(interiorborder_it);
509 g.interior[codim].finalize(interior_it);
510 g.send_overlapfront_overlapfront[codim].finalize(send_overlapfront_overlapfront_it,g.overlapfront[codim]);
511 g.recv_overlapfront_overlapfront[codim].finalize(recv_overlapfront_overlapfront_it,g.overlapfront[codim]);
512 g.send_overlap_overlapfront[codim].finalize(send_overlap_overlapfront_it,g.overlapfront[codim]);
513 g.recv_overlapfront_overlap[codim].finalize(recv_overlapfront_overlap_it,g.overlapfront[codim]);
514 g.send_interiorborder_interiorborder[codim].finalize(send_interiorborder_interiorborder_it,g.overlapfront[codim]);
515 g.recv_interiorborder_interiorborder[codim].finalize(recv_interiorborder_interiorborder_it,g.overlapfront[codim]);
516 g.send_interiorborder_overlapfront[codim].finalize(send_interiorborder_overlapfront_it,g.overlapfront[codim]);
517 g.recv_overlapfront_interiorborder[codim].finalize(recv_overlapfront_interiorborder_it,g.overlapfront[codim]);
530 struct mpifriendly_ygrid {
533 std::fill(origin.begin(), origin.end(), 0);
534 std::fill(size.begin(), size.end(), 0);
536 mpifriendly_ygrid (
const YGridComponent<Coordinates>& grid)
537 : origin(grid.origin()), size(grid.size())
553 std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
558 std::vector<YGridComponent<Coordinates> > send_recvgrid(_torus.
neighbors());
559 std::vector<YGridComponent<Coordinates> > recv_recvgrid(_torus.
neighbors());
560 std::vector<YGridComponent<Coordinates> > send_sendgrid(_torus.
neighbors());
561 std::vector<YGridComponent<Coordinates> > recv_sendgrid(_torus.
neighbors());
564 std::vector<mpifriendly_ygrid> mpifriendly_send_recvgrid(_torus.
neighbors());
565 std::vector<mpifriendly_ygrid> mpifriendly_recv_recvgrid(_torus.
neighbors());
566 std::vector<mpifriendly_ygrid> mpifriendly_send_sendgrid(_torus.
neighbors());
567 std::vector<mpifriendly_ygrid> mpifriendly_recv_sendgrid(_torus.
neighbors());
575 iTupel coord = _torus.
coord();
576 iTupel delta = i.delta();
578 for (
int k=0; k<dim; k++) nb[k] += delta[k];
580 std::fill(v.begin(), v.end(), 0);
582 for (
int k=0; k<dim; k++)
591 if (nb[k]>=_torus.
dims(k))
604 send_sendgrid[i.index()] = sendgrid.
move(v);
605 send_recvgrid[i.index()] = recvgrid.
move(v);
617 mpifriendly_send_sendgrid[i.index()] = mpifriendly_ygrid(send_sendgrid[i.index()]);
618 _torus.
send(i.rank(), &mpifriendly_send_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
623 _torus.
recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
631 mpifriendly_send_recvgrid[i.index()] = mpifriendly_ygrid(send_recvgrid[i.index()]);
632 _torus.
send(i.rank(), &mpifriendly_send_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
637 _torus.
recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
647 mpifriendly_ygrid yg = mpifriendly_recv_recvgrid[i.index()];
649 send_intersection.grid = sendgrid.
intersection(recv_recvgrid[i.index()]);
650 send_intersection.rank = i.rank();
651 send_intersection.distance = i.distance();
652 if (!send_intersection.grid.empty()) sendlist.push_front(send_intersection);
655 yg = mpifriendly_recv_sendgrid[i.index()];
657 recv_intersection.grid = recvgrid.
intersection(recv_sendgrid[i.index()]);
658 recv_intersection.rank = i.rank();
659 recv_intersection.distance = i.distance();
660 if(!recv_intersection.grid.empty()) recvlist.push_back(recv_intersection);
677 std::array<int, dim> sides;
679 for (
int i=0; i<dim; i++)
682 ((
begin()->overlap[0].dataBegin()->origin(i) == 0)+
683 (
begin()->overlap[0].dataBegin()->origin(i) +
begin()->overlap[0].dataBegin()->size(i)
688 for (
int k=0; k<dim; k++)
691 for (
int l=0; l<dim; l++)
694 offset *=
begin()->overlap[0].dataBegin()->size(l);
696 nBSegments += sides[k]*offset;
723 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
728 , leafIndexSet_(*this)
729 , _periodic(periodic)
738 for (std::size_t i=0; i<dim; i++)
739 _coarseSize[i] = coordinates.size(i);
742 _torus =
decltype(_torus)(
comm,tag,_coarseSize,lb);
745 std::fill(o.begin(), o.end(), 0);
746 iTupel o_interior(o);
747 iTupel s_interior(_coarseSize);
749 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
755 for (std::size_t i=0; i<dim; i++)
756 _L[i] = coordinates.size(i) * coordinates.meshsize(i,0);
761 for (
int i=0; i<dim; i++)
762 _L[i] = coordinates.coordinate(i,_coarseSize[i]) - coordinates.coordinate(i,0);
767 int mysteryFactor = (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value) ? 1 : 2;
770 for (
int i=0; i<dim; i++)
773 int toosmall = (s_interior[i] <= mysteryFactor * overlap) &&
774 (periodic[i] || (s_interior[i] != _coarseSize[i]));
777 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
779 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
780 " Note that this also holds for DOFs on subdomain boundaries."
781 " Increase grid elements or decrease overlap accordingly.");
788 iTupel s_overlap(s_interior);
789 for (
int i=0; i<dim; i++)
791 if ((o_interior[i] - overlap > 0) || (periodic[i]))
792 s_overlap[i] += overlap;
793 if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
794 s_overlap[i] += overlap;
797 FieldVector<ctype,dim> upperRightWithOverlap;
798 for (
int i=0; i<dim; i++)
799 upperRightWithOverlap[i] = coordinates.coordinate(i,0) + coordinates.meshsize(i,0) * s_overlap[i];
801 if constexpr (std::is_same_v<Coordinates,EquidistantCoordinates<ctype,dim>>)
807 this->
makelevel(coordinatesWithOverlap,periodic,o_interior,overlap);
810 if constexpr (std::is_same_v<Coordinates,EquidistantOffsetCoordinates<ctype,dim>>)
812 Dune::FieldVector<ctype,dim> lowerleft;
813 for (
int i=0; i<dim; i++)
814 lowerleft[i] = coordinates.origin(i);
820 this->
makelevel(coordinatesWithOverlap,periodic,o_interior,overlap);
824 if constexpr (std::is_same_v<Coordinates,TensorProductCoordinates<ctype,dim>>)
826 std::array<std::vector<ctype>,dim> newCoords;
827 std::array<int, dim> offset(o_interior);
830 for (
int i=0; i<dim; ++i)
833 std::size_t
begin = o_interior[i];
834 std::size_t
end =
begin + s_interior[i] + 1;
838 if (o_interior[i] - overlap > 0)
841 offset[i] -= overlap;
843 if (o_interior[i] + s_interior[i] + overlap < _coarseSize[i])
848 auto newCoordsIt = newCoords[i].begin();
849 for (std::size_t j=
begin; j<
end; j++)
851 *newCoordsIt = coordinates.coordinate(i, j);
857 if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
860 for (
int j=0; j<overlap; ++j)
861 newCoords[i].push_back(newCoords[i].back() - coordinates.coordinate(i,j) + coordinates.coordinate(i,j+1));
864 if ((periodic[i]) && (o_interior[i] - overlap <= 0))
866 offset[i] -= overlap;
869 std::size_t reverseCounter = coordinates.size(i);
870 for (
int j=0; j<overlap; ++j, --reverseCounter)
871 newCoords[i].insert(newCoords[i].
begin(), newCoords[i].front()
872 - coordinates.coordinate(i,reverseCounter) + coordinates.coordinate(i,reverseCounter-1));
879 this->
makelevel(coordinatesWithOverlap,periodic,o_interior,overlap);
893 template<
class C = Coordinates,
894 typename std::enable_if_t< std::is_same_v<C, EquidistantCoordinates<ctype,dim> >,
int> = 0>
896 std::array<
int, std::size_t{dim}> s,
897 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>{0ULL},
901 : ccobj(
comm), _torus(
comm,tag,s,lb), leafIndexSet_(*
this),
902 _L(L), _periodic(periodic), _coarseSize(s), _overlap(overlap),
903 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
908 std::fill(o.begin(), o.end(), 0);
909 iTupel o_interior(o);
910 iTupel s_interior(s);
916 for (
int i=0; i<dim; i++)
919 int toosmall = (s_interior[i] / 2 <=
overlap) &&
920 (periodic[i] || (s_interior[i] != s[i]));
923 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
925 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
926 " Note that this also holds for DOFs on subdomain boundaries."
927 " Increase grid elements or decrease overlap accordingly.");
931 iTupel s_overlap(s_interior);
932 for (
int i=0; i<dim; i++)
934 if ((o_interior[i] - overlap > 0) || (periodic[i]))
936 if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
940 FieldVector<ctype,dim> upperRightWithOverlap;
942 for (
int i=0; i<dim; i++)
943 upperRightWithOverlap[i] = (L[i] / s[i]) * s_overlap[i];
946 EquidistantCoordinates<ctype,dim> cc(upperRightWithOverlap,s_overlap);
949 makelevel(cc,periodic,o_interior,overlap);
963 template<
class C = Coordinates,
964 typename std::enable_if_t< std::is_same_v<C, EquidistantOffsetCoordinates<ctype,dim> >,
int> = 0>
966 Dune::FieldVector<ctype, dim> upperright,
967 std::array<
int, std::size_t{dim}> s,
968 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>(0ULL),
972 : ccobj(
comm), _torus(
comm,tag,s,lb), leafIndexSet_(*
this),
973 _L(upperright - lowerleft),
974 _periodic(periodic), _coarseSize(s), _overlap(overlap),
975 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
980 std::fill(o.begin(), o.end(), 0);
981 iTupel o_interior(o);
982 iTupel s_interior(s);
988 for (
int i=0; i<dim; i++)
991 int toosmall = (s_interior[i] / 2 <=
overlap) &&
992 (periodic[i] || (s_interior[i] != s[i]));
995 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
997 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
998 " Note that this also holds for DOFs on subdomain boundaries."
999 " Increase grid elements or decrease overlap accordingly.");
1003 iTupel s_overlap(s_interior);
1004 for (
int i=0; i<dim; i++)
1006 if ((o_interior[i] - overlap > 0) || (periodic[i]))
1008 if ((o_interior[i] + s_interior[i] + overlap <= _coarseSize[i]) || (periodic[i]))
1012 FieldVector<ctype,dim> upperRightWithOverlap;
1013 for (
int i=0; i<dim; i++)
1014 upperRightWithOverlap[i] = lowerleft[i]
1015 + s_overlap[i] * (upperright[i]-lowerleft[i]) / s[i];
1017 EquidistantOffsetCoordinates<ctype,dim> cc(lowerleft,upperRightWithOverlap,s_overlap);
1020 makelevel(cc,periodic,o_interior,overlap);
1032 template<
class C = Coordinates,
1033 typename std::enable_if_t< std::is_same_v<C, TensorProductCoordinates<ctype,dim> >,
int> = 0>
1034 YaspGrid (std::array<std::vector<ctype>, std::size_t{dim}> coords,
1035 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>(0ULL),
1039 : ccobj(
comm), _torus(
comm,tag,Dune::Yasp::sizeArray<dim>(coords),lb),
1040 leafIndexSet_(*
this), _periodic(periodic), _overlap(overlap),
1041 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
1044 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1049 for (
int i=0; i<dim; i++) {
1050 _coarseSize[i] = coords[i].size() - 1;
1051 _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
1055 std::fill(o.begin(), o.end(), 0);
1056 iTupel o_interior(o);
1057 iTupel s_interior(_coarseSize);
1059 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
1063 for (
int i=0; i<dim; i++)
1066 int toosmall = (s_interior[i] / 2 <=
overlap) &&
1067 (periodic[i] || (s_interior[i] != _coarseSize[i]));
1070 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
1072 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1073 " Note that this also holds for DOFs on subdomain boundaries."
1074 " Increase grid elements or decrease overlap accordingly.");
1079 std::array<std::vector<ctype>,dim> newcoords;
1080 std::array<int, dim> offset(o_interior);
1083 for (
int i=0; i<dim; ++i)
1086 typename std::vector<ctype>::iterator
begin = coords[i].begin() + o_interior[i];
1087 typename std::vector<ctype>::iterator
end =
begin + s_interior[i] + 1;
1091 if (o_interior[i] - overlap > 0)
1096 if (o_interior[i] + s_interior[i] + overlap < _coarseSize[i])
1105 if ((periodic[i]) && (o_interior[i] + s_interior[i] + overlap >= _coarseSize[i]))
1108 typename std::vector<ctype>::iterator it = coords[i].begin();
1110 newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
1113 if ((periodic[i]) && (o_interior[i] - overlap <= 0))
1118 typename std::vector<ctype>::iterator it = coords[i].end() - 1;
1120 newcoords[i].insert(newcoords[i].
begin(), newcoords[i].
front() - *it + *(--it));
1124 TensorProductCoordinates<ctype,dim> cc(newcoords, offset);
1127 makelevel(cc,periodic,o_interior,overlap);
1147 YaspGrid (std::array<std::vector<ctype>, std::size_t{dim}> coords,
1148 std::bitset<std::size_t{dim}> periodic,
1151 std::array<int,dim> coarseSize,
1153 : ccobj(
comm), _torus(
comm,tag,coarseSize,lb), leafIndexSet_(*
this),
1154 _periodic(periodic), _coarseSize(coarseSize), _overlap(overlap),
1155 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
1158 static_assert(std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >::value,
1159 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1162 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1164 for (
int i=0; i<dim; i++)
1165 _L[i] = coords[i][coords[i].
size() - 1] - coords[i][0];
1169 std::array<int,dim> o;
1170 std::fill(o.begin(), o.end(), 0);
1171 std::array<int,dim> o_interior(o);
1172 std::array<int,dim> s_interior(coarseSize);
1174 _torus.
partition(_torus.
rank(),o,coarseSize,o_interior,s_interior);
1177 std::array<int,dim> offset(o_interior);
1178 for (
int i=0; i<dim; i++)
1179 if ((periodic[i]) || (o_interior[i] > 0))
1180 offset[i] -= overlap;
1182 TensorProductCoordinates<ctype,dim> cc(coords, offset);
1185 makelevel(cc,periodic,o_interior,overlap);
1191 friend struct BackupRestoreFacility<
YaspGrid<dim,Coordinates> >;
1203 return _levels.size()-1;
1211 "Coarsening " << -refCount <<
" levels requested!");
1214 for (
int k=refCount; k<0; k++)
1218 _levels.back() = empty;
1222 indexsets.pop_back();
1226 for (
int k=0; k<refCount; k++)
1229 YGridLevel& cg = _levels[
maxLevel()];
1231 std::bitset<dim> ovlp_low(0ULL), ovlp_up(0ULL);
1232 for (
int i=0; i<dim; i++)
1234 if (cg.overlap[0].dataBegin()->origin(i) > 0 || _periodic[i])
1236 if (cg.overlap[0].dataBegin()->max(i) + 1 <
globalSize(i) || _periodic[i])
1240 Coordinates newcont(cg.coords.refine(ovlp_low, ovlp_up, cg.overlapSize, keep_ovlp));
1242 int overlap = (keep_ovlp) ? 2*cg.overlapSize : cg.overlapSize;
1246 for (
int i=0; i<dim; i++)
1247 o_interior[i] = 2*cg.interior[0].dataBegin()->origin(i);
1250 _levels.resize(_levels.size() + 1);
1251 makelevel(newcont,_periodic,o_interior,overlap);
1263 keep_ovlp = keepPhysicalOverlap;
1279 assert(adaptActive ==
false);
1280 if (e.level() !=
maxLevel())
return false;
1281 adaptRefCount = std::max(adaptRefCount, refCount);
1293 return ( e.level() ==
maxLevel() ) ? adaptRefCount : 0;
1300 return (adaptRefCount > 0);
1307 adaptRefCount =
comm().max(adaptRefCount);
1308 return (adaptRefCount < 0);
1314 adaptActive =
false;
1319 template<
int cd, PartitionIteratorType pitype>
1322 return levelbegin<cd,pitype>(level);
1326 template<
int cd, PartitionIteratorType pitype>
1329 return levelend<cd,pitype>(level);
1336 return levelbegin<cd,All_Partition>(level);
1343 return levelend<cd,All_Partition>(level);
1347 template<
int cd, PartitionIteratorType pitype>
1350 return levelbegin<cd,pitype>(
maxLevel());
1354 template<
int cd, PartitionIteratorType pitype>
1357 return levelend<cd,pitype>(
maxLevel());
1364 return levelbegin<cd,All_Partition>(
maxLevel());
1371 return levelend<cd,All_Partition>(
maxLevel());
1375 template <
typename Seed>
1376 typename Traits::template Codim<Seed::codimension>::Entity
1379 const int codim = Seed::codimension;
1386 return Entity(EntityImp(g,YIterator(g->overlapfront[codim],seed.impl().coord(),seed.impl().offset())));
1393 return g->overlapSize;
1400 return g->overlapSize;
1404 int ghostSize ([[maybe_unused]]
int level, [[maybe_unused]]
int codim)
const
1416 int size (
int level,
int codim)
const
1422 typedef typename std::array<YGridComponent<Coordinates>, StaticPower<2,dim>::power>::iterator DAI;
1423 for (DAI it = g->overlapfront[codim].dataBegin(); it != g->overlapfront[codim].dataEnd(); ++it)
1424 count += it->totalsize();
1436 int size (
int level, GeometryType type)
const
1438 return (type.isCube()) ?
size(level,dim-type.dim()) : 0;
1462 template<
class DataHandleImp,
class DataType>
1465 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,level);
1472 template<
class DataHandleImp,
class DataType>
1475 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,this->
maxLevel());
1482 template<
class DataHandle,
int codim>
1486 if (!data.contains(dim,codim))
return;
1489 typedef typename DataHandle::DataType DataType;
1500 sendlist = &g->send_interiorborder_interiorborder[codim];
1501 recvlist = &g->recv_interiorborder_interiorborder[codim];
1505 sendlist = &g->send_interiorborder_overlapfront[codim];
1506 recvlist = &g->recv_overlapfront_interiorborder[codim];
1510 sendlist = &g->send_overlap_overlapfront[codim];
1511 recvlist = &g->recv_overlapfront_overlap[codim];
1515 sendlist = &g->send_overlapfront_overlapfront[codim];
1516 recvlist = &g->recv_overlapfront_overlapfront[codim];
1526 std::vector<int> send_size(sendlist->
size(),-1);
1527 std::vector<int> recv_size(recvlist->
size(),-1);
1528 std::vector<size_t*> send_sizes(sendlist->
size(),
static_cast<size_t*
>(0));
1529 std::vector<size_t*> recv_sizes(recvlist->
size(),
static_cast<size_t*
>(0));
1534 if (data.fixedSize(dim,codim))
1538 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1542 send_size[cnt] = is->grid.totalsize() * data.size(*it);
1546 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1550 recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1558 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1561 size_t *buf =
new size_t[is->grid.totalsize()];
1562 send_sizes[cnt] = buf;
1565 int i=0;
size_t n=0;
1570 for ( ; it!=itend; ++it)
1572 buf[i] = data.size(*it);
1581 torus().
send(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1587 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1590 size_t *buf =
new size_t[is->grid.totalsize()];
1591 recv_sizes[cnt] = buf;
1594 torus().
recv(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1603 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1605 delete[] send_sizes[cnt];
1606 send_sizes[cnt] = 0;
1612 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1615 size_t *buf = recv_sizes[cnt];
1619 for (
int i=0; i<is->grid.totalsize(); ++i)
1630 std::vector<DataType*> sends(sendlist->
size(),
static_cast<DataType*
>(0));
1632 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1635 DataType *buf =
new DataType[send_size[cnt]];
1641 MessageBuffer<DataType> mb(buf);
1648 for ( ; it!=itend; ++it)
1649 data.gather(mb,*it);
1652 torus().
send(is->rank,buf,send_size[cnt]*
sizeof(DataType));
1657 std::vector<DataType*> recvs(recvlist->
size(),
static_cast<DataType*
>(0));
1659 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1662 DataType *buf =
new DataType[recv_size[cnt]];
1668 torus().
recv(is->rank,buf,recv_size[cnt]*
sizeof(DataType));
1677 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1679 delete[] sends[cnt];
1686 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1689 DataType *buf = recvs[cnt];
1692 MessageBuffer<DataType> mb(buf);
1695 if (data.fixedSize(dim,codim))
1699 size_t n=data.size(*it);
1702 for ( ; it!=itend; ++it)
1703 data.scatter(mb,*it,n);
1708 size_t *sbuf = recv_sizes[cnt];
1713 for ( ; it!=itend; ++it)
1714 data.scatter(mb,*it,sbuf[i++]);
1727 return theglobalidset;
1732 return theglobalidset;
1737 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1738 return *(indexsets[level]);
1743 return leafIndexSet_;
1766 friend class
Dune::YaspEntity<0, dim, const Dune::YaspGrid<dim, Coordinates> >;
1768 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1772 class MessageBuffer {
1775 MessageBuffer (DT *p)
1784 void write (
const Y& data)
1786 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1792 void read (Y& data)
const
1794 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1805 template<
int cd, PartitionIteratorType pitype>
1809 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1820 return levelend <cd, pitype> (level);
1822 DUNE_THROW(
GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1826 template<
int cd, PartitionIteratorType pitype>
1827 YaspLevelIterator<cd,pitype,GridImp> levelend (
int level)
const
1830 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1833 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].end());
1835 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].end());
1837 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].end());
1839 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].end());
1841 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1846 Torus<CollectiveCommunicationType,dim> _torus;
1848 std::vector< std::shared_ptr< YaspIndexSet<const YaspGrid<dim,Coordinates>,
false > > > indexsets;
1849 YaspIndexSet<const YaspGrid<dim,Coordinates>,
true> leafIndexSet_;
1850 YaspGlobalIdSet<const YaspGrid<dim,Coordinates> > theglobalidset;
1852 Dune::FieldVector<ctype, dim> _L;
1854 std::bitset<dim> _periodic;
1856 ReservedVector<YGridLevel,32> _levels;
1863#if __cpp_deduction_guides >= 201611
1865 template<
typename ctype,
int dim>
1866 YaspGrid(FieldVector<ctype, dim>,
1867 std::array<
int, std::size_t{dim}>,
1868 std::bitset<std::size_t{dim}> = std::bitset<std::size_t{dim}>{0ULL},
1871 const YLoadBalance<dim>* = YaspGrid< dim, EquidistantCoordinates<ctype, dim> >::defaultLoadbalancer())
1872 -> YaspGrid< dim, EquidistantCoordinates<ctype, dim> >;
1874 template<
typename ctype,
int dim>
1875 YaspGrid(FieldVector<ctype, dim>,
1876 FieldVector<ctype, dim>,
1877 std::array<
int, std::size_t{dim}>,
1878 std::bitset<std::size_t{dim}> = std::bitset<std::size_t{dim}>{0ULL},
1881 const YLoadBalance<dim>* = YaspGrid< dim, EquidistantOffsetCoordinates<ctype, dim> >::defaultLoadbalancer())
1882 -> YaspGrid< dim, EquidistantOffsetCoordinates<ctype, dim> >;
1884 template<
typename ctype, std::
size_t dim>
1885 YaspGrid(std::array<std::vector<ctype>, dim>,
1886 std::bitset<dim> = std::bitset<dim>{0ULL},
1889 const YLoadBalance<
int{dim}>* = YaspGrid<
int{dim}, TensorProductCoordinates<ctype,
int{dim}> >::defaultLoadbalancer())
1890 -> YaspGrid<
int{dim}, TensorProductCoordinates<ctype,
int{dim}> >;
1894 template <
int d,
class CC>
1899 s <<
"[" << rank <<
"]:" <<
" YaspGrid maxlevel=" << grid.
maxLevel() << std::endl;
1901 s <<
"Printing the torus: " <<std::endl;
1902 s << grid.
torus() << std::endl;
1906 s <<
"[" << rank <<
"]: " << std::endl;
1907 s <<
"[" << rank <<
"]: " <<
"==========================================" << std::endl;
1908 s <<
"[" << rank <<
"]: " <<
"level=" << g->level() << std::endl;
1910 for (
int codim = 0; codim < d + 1; ++codim)
1912 s <<
"[" << rank <<
"]: " <<
"overlapfront[" << codim <<
"]: " << g->overlapfront[codim] << std::endl;
1913 s <<
"[" << rank <<
"]: " <<
"overlap[" << codim <<
"]: " << g->overlap[codim] << std::endl;
1914 s <<
"[" << rank <<
"]: " <<
"interiorborder[" << codim <<
"]: " << g->interiorborder[codim] << std::endl;
1915 s <<
"[" << rank <<
"]: " <<
"interior[" << codim <<
"]: " << g->interior[codim] << std::endl;
1918 for (I i=g->send_overlapfront_overlapfront[codim].begin();
1919 i!=g->send_overlapfront_overlapfront[codim].end(); ++i)
1920 s <<
"[" << rank <<
"]: " <<
" s_of_of[" << codim <<
"] to rank "
1921 << i->rank <<
" " << i->grid << std::endl;
1923 for (I i=g->recv_overlapfront_overlapfront[codim].begin();
1924 i!=g->recv_overlapfront_overlapfront[codim].end(); ++i)
1925 s <<
"[" << rank <<
"]: " <<
" r_of_of[" << codim <<
"] to rank "
1926 << i->rank <<
" " << i->grid << std::endl;
1928 for (I i=g->send_overlap_overlapfront[codim].begin();
1929 i!=g->send_overlap_overlapfront[codim].end(); ++i)
1930 s <<
"[" << rank <<
"]: " <<
" s_o_of[" << codim <<
"] to rank "
1931 << i->rank <<
" " << i->grid << std::endl;
1933 for (I i=g->recv_overlapfront_overlap[codim].begin();
1934 i!=g->recv_overlapfront_overlap[codim].end(); ++i)
1935 s <<
"[" << rank <<
"]: " <<
" r_of_o[" << codim <<
"] to rank "
1936 << i->rank <<
" " << i->grid << std::endl;
1938 for (I i=g->send_interiorborder_interiorborder[codim].begin();
1939 i!=g->send_interiorborder_interiorborder[codim].end(); ++i)
1940 s <<
"[" << rank <<
"]: " <<
" s_ib_ib[" << codim <<
"] to rank "
1941 << i->rank <<
" " << i->grid << std::endl;
1943 for (I i=g->recv_interiorborder_interiorborder[codim].begin();
1944 i!=g->recv_interiorborder_interiorborder[codim].end(); ++i)
1945 s <<
"[" << rank <<
"]: " <<
" r_ib_ib[" << codim <<
"] to rank "
1946 << i->rank <<
" " << i->grid << std::endl;
1948 for (I i=g->send_interiorborder_overlapfront[codim].begin();
1949 i!=g->send_interiorborder_overlapfront[codim].end(); ++i)
1950 s <<
"[" << rank <<
"]: " <<
" s_ib_of[" << codim <<
"] to rank "
1951 << i->rank <<
" " << i->grid << std::endl;
1953 for (I i=g->recv_overlapfront_interiorborder[codim].begin();
1954 i!=g->recv_overlapfront_interiorborder[codim].end(); ++i)
1955 s <<
"[" << rank <<
"]: " <<
" r_of_ib[" << codim <<
"] to rank "
1956 << i->rank <<
" " << i->grid << std::endl;
1965 namespace Capabilities
1975 template<
int dim,
class Coordinates>
1978 static const bool v =
true;
1984 template<
int dim,
class Coordinates>
1987 static const bool v =
true;
1988 static const unsigned int topologyId = GeometryTypes::cube(dim).id();
1994 template<
int dim,
class Coordinates>
1997 static const bool v =
true;
2003 template<
int dim,
class Coordinates,
int codim>
2006 static const bool v =
true;
2013 template<
int dim,
class Coordinates,
int codim>
2016 static const bool v =
true;
2022 template<
int dim,
int codim,
class Coordinates>
2025 static const bool v =
true;
2031 template<
int dim,
class Coordinates>
2034 static const bool v =
true;
2040 template<
int dim,
class Coordinates>
2043 static const bool v =
true;
Specialization of the StructuredGridFactory class for YaspGrid.
The YaspLevelIterator class.
This provides a YGrid, the elemental component of the yaspgrid implementation.
The YaspIntersectionIterator class.
The YaspIntersection class.
level-wise, non-persistent, consecutive indices for YaspGrid
This file provides the infrastructure for toroidal communication in YaspGrid.
The YaspEntitySeed class.
This provides container classes for the coordinates to be used in YaspGrid Upon implementation of the...
the YaspEntity class and its specializations
Specialization of the PersistentContainer for YaspGrid.
The YaspGeometry class and its specializations.
Describes the parallel communication interface class for MessageBuffers and DataHandles.
Provides base classes for index and id sets.
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
void swap(Dune::PersistentContainer< G, T > &a, Dune::PersistentContainer< G, T > &b)
Definition: utility/persistentcontainer.hh:81
Include standard header files.
Definition: agrid.hh:58
const int yaspgrid_level_bits
Definition: yaspgrid.hh:45
CollectiveCommunication< MPI_Comm > YaspCollectiveCommunication
Definition: yaspgrid.hh:82
const int yaspgrid_dim_bits
Definition: yaspgrid.hh:44
constexpr Overlap overlap
PartitionSet for the overlap partition.
Definition: partitionset.hh:276
constexpr Front front
PartitionSet for the front partition.
Definition: partitionset.hh:279
bool checkIfMonotonous(const std::array< std::vector< ctype >, dim > &coords)
Definition: coordinates.hh:370
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:162
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=hasEntity<codim>:...
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:851
Index Set Interface base class.
Definition: indexidset.hh:76
Id Set Interface.
Definition: indexidset.hh:450
A Traits struct that collects all associated types of one implementation.
Definition: common/grid.hh:414
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:984
[ provides Dune::Grid ]
Definition: yaspgrid.hh:160
YaspGridFamily< dim, Coordinates > GridFamily
the GridFamily of this grid
Definition: yaspgrid.hh:706
YaspIndexSet< YaspGrid< dim, Coordinates >, false > LevelIndexSetType
Definition: yaspgrid.hh:711
int ghostSize(int level, int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1404
YaspGrid(Dune::FieldVector< ctype, dim > L, std::array< int, std::size_t{dim}> s, std::bitset< std::size_t{dim}> periodic=std::bitset< std::size_t{dim}>{0ULL}, int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:895
const Traits::LeafIndexSet & leafIndexSet() const
Definition: yaspgrid.hh:1741
void init()
Definition: yaspgrid.hh:668
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lbegin(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1334
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: yaspgrid.hh:1436
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Definition: yaspgrid.hh:1377
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1362
const Traits::GlobalIdSet & globalIdSet() const
Definition: yaspgrid.hh:1725
const Traits::LocalIdSet & localIdSet() const
Definition: yaspgrid.hh:1730
YaspCollectiveCommunication CollectiveCommunicationType
Definition: yaspgrid.hh:173
YaspGridFamily< dim, Coordinates >::Traits Traits
Definition: yaspgrid.hh:708
bool getRefineOption() const
Definition: yaspgrid.hh:276
void globalRefine(int refCount)
refine the grid refCount times.
Definition: yaspgrid.hh:1207
const Torus< CollectiveCommunicationType, dim > & torus() const
return reference to torus
Definition: yaspgrid.hh:238
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
returns adaptation mark for given entity
Definition: yaspgrid.hh:1291
void boundarysegmentssize()
Definition: yaspgrid.hh:674
YaspIndexSet< YaspGrid< dim, Coordinates >, true > LeafIndexSetType
Definition: yaspgrid.hh:712
int ghostSize(int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1410
int globalSize(int i) const
return number of cells on finest level in given direction on all processors
Definition: yaspgrid.hh:244
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1355
void postAdapt()
clean up some markers
Definition: yaspgrid.hh:1312
const Dune::FieldVector< ctype, dim > & domainSize() const
returns the size of the physical domain
Definition: yaspgrid.hh:1454
YGridLevelIterator end() const
return iterator pointing to one past the finest level
Definition: yaspgrid.hh:299
friend class Entity
Definition: yaspgrid.hh:1769
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: yaspgrid.hh:1442
int maxLevel() const
Definition: yaspgrid.hh:1201
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1369
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1463
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:552
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1348
bool preAdapt()
returns true, if the grid will be coarsened
Definition: yaspgrid.hh:1304
iTupel levelSize(int l) const
return size vector of the grid (in cells) on level l
Definition: yaspgrid.hh:262
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:1277
int overlapSize(int level, int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1390
const CollectiveCommunicationType & comm() const
return a collective communication object
Definition: yaspgrid.hh:1748
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir) const
Definition: yaspgrid.hh:1473
int size(int codim) const
number of leaf entities per codim in this process
Definition: yaspgrid.hh:1430
bool isPeriodic(int i) const
return whether the grid is periodic in direction i
Definition: yaspgrid.hh:271
const Traits::LevelIndexSet & levelIndexSet(int level) const
Definition: yaspgrid.hh:1735
void refineOptions(bool keepPhysicalOverlap)
set options for refinement
Definition: yaspgrid.hh:1261
int size(int level, int codim) const
number of entities per level and codim in this process
Definition: yaspgrid.hh:1416
YaspGlobalIdSet< YaspGrid< dim, Coordinates > > GlobalIdSetType
Definition: yaspgrid.hh:713
bool adapt()
map adapt to global refine
Definition: yaspgrid.hh:1297
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lend(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1341
ReservedVector< YGridLevel, 32 >::const_iterator YGridLevelIterator
Iterator over the grid levels.
Definition: yaspgrid.hh:282
size_t numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: yaspgrid.hh:1448
void communicateCodim(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1483
void makelevel(const Coordinates &coords, std::bitset< dim > periodic, iTupel o_interior, int overlap)
Make a new YGridLevel structure.
Definition: yaspgrid.hh:319
bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim > PersistentIndexType
Definition: yaspgrid.hh:703
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:1327
static const YLoadBalanceDefault< dim > * defaultLoadbalancer()
Definition: yaspgrid.hh:305
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lbegin(int level) const
one past the end on this level
Definition: yaspgrid.hh:1320
YaspGrid(std::array< std::vector< ctype >, std::size_t{dim}> coords, std::bitset< std::size_t{dim}> periodic=std::bitset< std::size_t{dim}>(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Standard constructor for a tensorproduct YaspGrid.
Definition: yaspgrid.hh:1034
YaspGrid(Dune::FieldVector< ctype, dim > lowerleft, Dune::FieldVector< ctype, dim > upperright, std::array< int, std::size_t{dim}> s, std::bitset< std::size_t{dim}> periodic=std::bitset< std::size_t{dim}>(0ULL), int overlap=1, CollectiveCommunicationType comm=CollectiveCommunicationType(), const YLoadBalance< dim > *lb=defaultLoadbalancer())
Definition: yaspgrid.hh:965
iTupel globalSize() const
return number of cells on finest level on all processors
Definition: yaspgrid.hh:250
int overlapSize(int odim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1397
const YaspGrid< dim, Coordinates > GridImp
Definition: yaspgrid.hh:666
YGridLevelIterator begin(int i) const
return iterator pointing to given level
Definition: yaspgrid.hh:291
Coordinates::ctype ctype
Type used for coordinates.
Definition: yaspgrid.hh:172
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:722
int levelSize(int l, int i) const
return size of the grid (in cells) on level l in direction i
Definition: yaspgrid.hh:256
YGridLevelIterator begin() const
return iterator pointing to coarsest level
Definition: yaspgrid.hh:285
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:89
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:112
YaspCollectiveCommunication CCType
Definition: yaspgrid.hh:90
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
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
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
const iTupel & dims() const
return dimensions of torus
Definition: torus.hh:110
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.