5#ifndef DUNE_GRID_YASPGRID_HH
6#define DUNE_GRID_YASPGRID_HH
19#include <dune/common/hybridutilities.hh>
20#include <dune/common/bigunsignedint.hh>
21#include <dune/common/math.hh>
22#include <dune/common/typetraits.hh>
23#include <dune/common/reservedvector.hh>
24#include <dune/common/parallel/communication.hh>
25#include <dune/common/parallel/mpihelper.hh>
26#include <dune/geometry/axisalignedcubegeometry.hh>
27#include <dune/geometry/type.hh>
33#include <dune/common/parallel/mpicommunication.hh>
54 template<
int dim,
class Coordinates>
class YaspGrid;
55 template<
int mydim,
int cdim,
class Gr
idImp>
class YaspGeometry;
56 template<
int codim,
int dim,
class Gr
idImp>
class YaspEntity;
58 template<
int codim, PartitionIteratorType pitype,
class Gr
idImp>
class YaspLevelIterator;
62 template<
class Gr
idImp,
bool isLeafIndexSet>
class YaspIndexSet;
90 template<
int dim,
class Coordinates>
109 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
111 bigunsignedint<dim*yaspgrid_dim_bits+yaspgrid_level_bits+dim>,
119 template<
int dim,
int codim>
120 struct YaspCommunicateMeta {
121 template<
class G,
class DataHandle>
124 if (data.contains(dim,codim))
126 g.template communicateCodim<DataHandle,codim>(data,iftype,dir,level);
128 YaspCommunicateMeta<dim,codim-1>::comm(g,data,iftype,dir,level);
133 struct YaspCommunicateMeta<dim,0> {
134 template<
class G,
class DataHandle>
137 if (data.contains(dim,0))
138 g.template communicateCodim<DataHandle,0>(data,iftype,dir,level);
160 template<
int dim,
class Coordinates = Equ
idistantCoordinates<
double, dim> >
165 template<
int, PartitionIteratorType,
typename>
176 typedef typename Coordinates::ctype
ctype;
197 std::array<YGrid, dim+1> overlapfront;
198 std::array<YGridComponent<Coordinates>, Dune::power(2,dim)> overlapfront_data;
199 std::array<YGrid, dim+1>
overlap;
200 std::array<YGridComponent<Coordinates>, Dune::power(2,dim)> overlap_data;
201 std::array<YGrid, dim+1> interiorborder;
202 std::array<YGridComponent<Coordinates>, Dune::power(2,dim)> interiorborder_data;
204 std::array<YGridComponent<Coordinates>, Dune::power(2,dim)> interior_data;
206 std::array<YGridList<Coordinates>,dim+1> send_overlapfront_overlapfront;
207 std::array<std::deque<Intersection>, Dune::power(2,dim)> send_overlapfront_overlapfront_data;
208 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlapfront;
209 std::array<std::deque<Intersection>, Dune::power(2,dim)> recv_overlapfront_overlapfront_data;
211 std::array<YGridList<Coordinates>,dim+1> send_overlap_overlapfront;
212 std::array<std::deque<Intersection>, Dune::power(2,dim)> send_overlap_overlapfront_data;
213 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_overlap;
214 std::array<std::deque<Intersection>, Dune::power(2,dim)> recv_overlapfront_overlap_data;
216 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_interiorborder;
217 std::array<std::deque<Intersection>, Dune::power(2,dim)> send_interiorborder_interiorborder_data;
218 std::array<YGridList<Coordinates>,dim+1> recv_interiorborder_interiorborder;
219 std::array<std::deque<Intersection>, Dune::power(2,dim)> recv_interiorborder_interiorborder_data;
221 std::array<YGridList<Coordinates>,dim+1> send_interiorborder_overlapfront;
222 std::array<std::deque<Intersection>, Dune::power(2,dim)> send_interiorborder_overlapfront_data;
223 std::array<YGridList<Coordinates>,dim+1> recv_overlapfront_interiorborder;
224 std::array<std::deque<Intersection>, Dune::power(2,dim)> recv_overlapfront_interiorborder_data;
236 typedef std::array<int, dim> iTupel;
237 typedef FieldVector<ctype, dim> fTupel;
264 return _coarseSize[i] * (1 << l);
271 for (
int i=0; i<dim; ++i)
293 return _levels.begin();
300 DUNE_THROW(
GridError,
"level not existing");
301 return std::next(_levels.begin(), i);
307 return _levels.end();
311 [[deprecated(
"use defaultPartitioner")]]
333 void makelevel (
const Coordinates& coords, std::bitset<dim> periodic, iTupel o_interior,
int overlap)
335 YGridLevel& g = _levels.back();
340 g.keepOverlap = keep_ovlp;
345 typename std::array<YGridComponent<Coordinates>, power(2,dim)>::iterator overlapfront_it = g.overlapfront_data.begin();
346 typename std::array<YGridComponent<Coordinates>, power(2,dim)>::iterator overlap_it = g.overlap_data.begin();
347 typename std::array<YGridComponent<Coordinates>, power(2,dim)>::iterator interiorborder_it = g.interiorborder_data.begin();
348 typename std::array<YGridComponent<Coordinates>, power(2,dim)>::iterator interior_it = g.interior_data.begin();
350 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
351 send_overlapfront_overlapfront_it = g.send_overlapfront_overlapfront_data.begin();
352 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
353 recv_overlapfront_overlapfront_it = g.recv_overlapfront_overlapfront_data.begin();
355 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
356 send_overlap_overlapfront_it = g.send_overlap_overlapfront_data.begin();
357 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
358 recv_overlapfront_overlap_it = g.recv_overlapfront_overlap_data.begin();
360 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
361 send_interiorborder_interiorborder_it = g.send_interiorborder_interiorborder_data.begin();
362 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
363 recv_interiorborder_interiorborder_it = g.recv_interiorborder_interiorborder_data.begin();
365 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
366 send_interiorborder_overlapfront_it = g.send_interiorborder_overlapfront_data.begin();
367 typename std::array<std::deque<Intersection>, power(2,dim)>::iterator
368 recv_overlapfront_interiorborder_it = g.recv_overlapfront_interiorborder_data.begin();
371 std::array<int,dim> n;
372 std::fill(n.begin(), n.end(), 0);
375 std::bitset<dim> ovlp_low(0ULL);
376 std::bitset<dim> ovlp_up(0ULL);
382 for (
int i=0; i<dim; i++)
386 s_overlap[i] = g.coords.size(i);
391 o_overlap[i] = o_interior[i]-
overlap;
398 if (o_interior[i] -
overlap < 0)
402 o_overlap[i] = o_interior[i] -
overlap;
407 if (o_overlap[i] + g.coords.size(i) <
globalSize(i))
412 for (
unsigned int codim = 0; codim < dim + 1; codim++)
415 g.overlapfront[codim].setBegin(overlapfront_it);
416 g.overlap[codim].setBegin(overlap_it);
417 g.interiorborder[codim].setBegin(interiorborder_it);
418 g.interior[codim].setBegin(interior_it);
419 g.send_overlapfront_overlapfront[codim].setBegin(send_overlapfront_overlapfront_it);
420 g.recv_overlapfront_overlapfront[codim].setBegin(recv_overlapfront_overlapfront_it);
421 g.send_overlap_overlapfront[codim].setBegin(send_overlap_overlapfront_it);
422 g.recv_overlapfront_overlap[codim].setBegin(recv_overlapfront_overlap_it);
423 g.send_interiorborder_interiorborder[codim].setBegin(send_interiorborder_interiorborder_it);
424 g.recv_interiorborder_interiorborder[codim].setBegin(recv_interiorborder_interiorborder_it);
425 g.send_interiorborder_overlapfront[codim].setBegin(send_interiorborder_overlapfront_it);
426 g.recv_overlapfront_interiorborder[codim].setBegin(recv_overlapfront_interiorborder_it);
429 for (
unsigned int index = 0; index < (1<<dim); index++)
432 std::bitset<dim> r(index);
433 if (r.count() != dim-codim)
437 std::array<int,dim> origin(o_overlap);
438 std::array<int,dim>
size(s_overlap);
442 for (
int i=0; i<dim; i++)
448 for (
int i=0; i<dim; i++)
464 for (
int i=0; i<dim; i++)
486 for (
int i=0; i<dim; i++)
501 intersections(*overlapfront_it,*overlapfront_it,*send_overlapfront_overlapfront_it, *recv_overlapfront_overlapfront_it);
502 intersections(*overlap_it,*overlapfront_it,*send_overlap_overlapfront_it, *recv_overlapfront_overlap_it);
503 intersections(*interiorborder_it,*interiorborder_it,*send_interiorborder_interiorborder_it,*recv_interiorborder_interiorborder_it);
504 intersections(*interiorborder_it,*overlapfront_it,*send_interiorborder_overlapfront_it,*recv_overlapfront_interiorborder_it);
511 ++send_overlapfront_overlapfront_it;
512 ++recv_overlapfront_overlapfront_it;
513 ++send_overlap_overlapfront_it;
514 ++recv_overlapfront_overlap_it;
515 ++send_interiorborder_interiorborder_it;
516 ++recv_interiorborder_interiorborder_it;
517 ++send_interiorborder_overlapfront_it;
518 ++recv_overlapfront_interiorborder_it;
522 g.overlapfront[codim].finalize(overlapfront_it);
523 g.overlap[codim].finalize(overlap_it);
524 g.interiorborder[codim].finalize(interiorborder_it);
525 g.interior[codim].finalize(interior_it);
526 g.send_overlapfront_overlapfront[codim].finalize(send_overlapfront_overlapfront_it,g.overlapfront[codim]);
527 g.recv_overlapfront_overlapfront[codim].finalize(recv_overlapfront_overlapfront_it,g.overlapfront[codim]);
528 g.send_overlap_overlapfront[codim].finalize(send_overlap_overlapfront_it,g.overlapfront[codim]);
529 g.recv_overlapfront_overlap[codim].finalize(recv_overlapfront_overlap_it,g.overlapfront[codim]);
530 g.send_interiorborder_interiorborder[codim].finalize(send_interiorborder_interiorborder_it,g.overlapfront[codim]);
531 g.recv_interiorborder_interiorborder[codim].finalize(recv_interiorborder_interiorborder_it,g.overlapfront[codim]);
532 g.send_interiorborder_overlapfront[codim].finalize(send_interiorborder_overlapfront_it,g.overlapfront[codim]);
533 g.recv_overlapfront_interiorborder[codim].finalize(recv_overlapfront_interiorborder_it,g.overlapfront[codim]);
546 struct mpifriendly_ygrid {
549 std::fill(origin.begin(), origin.end(), 0);
550 std::fill(
size.begin(),
size.end(), 0);
552 mpifriendly_ygrid (
const YGridComponent<Coordinates>& grid)
553 : origin(grid.origin()),
size(grid.
size())
569 std::deque<Intersection>& sendlist, std::deque<Intersection>& recvlist)
574 std::vector<YGridComponent<Coordinates> > send_recvgrid(_torus.
neighbors());
575 std::vector<YGridComponent<Coordinates> > recv_recvgrid(_torus.
neighbors());
576 std::vector<YGridComponent<Coordinates> > send_sendgrid(_torus.
neighbors());
577 std::vector<YGridComponent<Coordinates> > recv_sendgrid(_torus.
neighbors());
580 std::vector<mpifriendly_ygrid> mpifriendly_send_recvgrid(_torus.
neighbors());
581 std::vector<mpifriendly_ygrid> mpifriendly_recv_recvgrid(_torus.
neighbors());
582 std::vector<mpifriendly_ygrid> mpifriendly_send_sendgrid(_torus.
neighbors());
583 std::vector<mpifriendly_ygrid> mpifriendly_recv_sendgrid(_torus.
neighbors());
591 iTupel coord = _torus.
coord();
592 iTupel delta = i.delta();
594 for (
int k=0; k<dim; k++) nb[k] += delta[k];
596 std::fill(v.begin(), v.end(), 0);
598 for (
int k=0; k<dim; k++)
607 if (nb[k]>=_torus.
dims(k))
620 send_sendgrid[i.index()] = sendgrid.
move(v);
621 send_recvgrid[i.index()] = recvgrid.
move(v);
633 mpifriendly_send_sendgrid[i.index()] = mpifriendly_ygrid(send_sendgrid[i.index()]);
634 _torus.
send(i.rank(), &mpifriendly_send_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
639 _torus.
recv(i.rank(), &mpifriendly_recv_sendgrid[i.index()],
sizeof(mpifriendly_ygrid));
647 mpifriendly_send_recvgrid[i.index()] = mpifriendly_ygrid(send_recvgrid[i.index()]);
648 _torus.
send(i.rank(), &mpifriendly_send_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
653 _torus.
recv(i.rank(), &mpifriendly_recv_recvgrid[i.index()],
sizeof(mpifriendly_ygrid));
663 mpifriendly_ygrid yg = mpifriendly_recv_recvgrid[i.index()];
665 send_intersection.grid = sendgrid.
intersection(recv_recvgrid[i.index()]);
666 send_intersection.rank = i.rank();
667 send_intersection.distance = i.distance();
668 if (!send_intersection.grid.empty()) sendlist.push_front(send_intersection);
671 yg = mpifriendly_recv_sendgrid[i.index()];
673 recv_intersection.grid = recvgrid.
intersection(recv_sendgrid[i.index()]);
674 recv_intersection.rank = i.rank();
675 recv_intersection.distance = i.distance();
676 if(!recv_intersection.grid.empty()) recvlist.push_back(recv_intersection);
693 std::array<int, dim> sides;
695 for (
int i=0; i<dim; i++)
698 ((
begin()->overlap[0].dataBegin()->origin(i) == 0)+
699 (
begin()->overlap[0].dataBegin()->origin(i) +
begin()->overlap[0].dataBegin()->size(i)
704 for (
int k=0; k<dim; k++)
707 for (
int l=0; l<dim; l++)
710 offset *=
begin()->overlap[0].dataBegin()->size(l);
712 nBSegments += sides[k]*offset;
739 std::bitset<dim> periodic = std::bitset<dim>(0ULL),
744 , leafIndexSet_(*this)
745 , _periodic(periodic)
754 for (std::size_t i=0; i<dim; i++)
755 _coarseSize[i] = coordinates.size(i);
758 _torus =
decltype(_torus)(
comm,tag,_coarseSize,
overlap,partitioner);
761 std::fill(o.begin(), o.end(), 0);
762 iTupel o_interior(o);
763 iTupel s_interior(_coarseSize);
765 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
771 for (std::size_t i=0; i<dim; i++)
772 _L[i] = coordinates.size(i) * coordinates.meshsize(i,0);
777 for (
int i=0; i<dim; i++)
778 _L[i] = coordinates.coordinate(i,_coarseSize[i]) - coordinates.coordinate(i,0);
783 int mysteryFactor = (std::is_same<Coordinates,EquidistantCoordinates<ctype,dim> >::value) ? 1 : 2;
786 for (
int i=0; i<dim; i++)
789 int toosmall = (s_interior[i] <= mysteryFactor *
overlap) &&
790 (periodic[i] || (s_interior[i] != _coarseSize[i]));
793 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
795 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
796 " Note that this also holds for DOFs on subdomain boundaries."
797 " Increase grid elements or decrease overlap accordingly.");
804 iTupel s_overlap(s_interior);
805 for (
int i=0; i<dim; i++)
807 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
809 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
813 FieldVector<ctype,dim> upperRightWithOverlap;
814 for (
int i=0; i<dim; i++)
815 upperRightWithOverlap[i] = coordinates.coordinate(i,0) + coordinates.meshsize(i,0) * s_overlap[i];
817 if constexpr (std::is_same_v<Coordinates,EquidistantCoordinates<ctype,dim>>)
826 if constexpr (std::is_same_v<Coordinates,EquidistantOffsetCoordinates<ctype,dim>>)
828 Dune::FieldVector<ctype,dim> lowerleft;
829 for (
int i=0; i<dim; i++)
830 lowerleft[i] = coordinates.origin(i);
840 if constexpr (std::is_same_v<Coordinates,TensorProductCoordinates<ctype,dim>>)
842 std::array<std::vector<ctype>,dim> newCoords;
843 std::array<int, dim> offset(o_interior);
846 for (
int i=0; i<dim; ++i)
849 std::size_t
begin = o_interior[i];
850 std::size_t
end =
begin + s_interior[i] + 1;
854 if (o_interior[i] -
overlap > 0)
859 if (o_interior[i] + s_interior[i] +
overlap < _coarseSize[i])
864 auto newCoordsIt = newCoords[i].begin();
865 for (std::size_t j=
begin; j<
end; j++)
867 *newCoordsIt = coordinates.coordinate(i, j);
873 if ((periodic[i]) && (o_interior[i] + s_interior[i] +
overlap >= _coarseSize[i]))
877 newCoords[i].push_back(newCoords[i].back() - coordinates.coordinate(i,j) + coordinates.coordinate(i,j+1));
880 if ((periodic[i]) && (o_interior[i] -
overlap <= 0))
885 std::size_t reverseCounter = coordinates.size(i);
886 for (
int j=0; j<
overlap; ++j, --reverseCounter)
887 newCoords[i].insert(newCoords[i].
begin(), newCoords[i].
front()
888 - coordinates.coordinate(i,reverseCounter) + coordinates.coordinate(i,reverseCounter-1));
909 template<
class C = Coordinates,
910 typename std::enable_if_t< std::is_same_v<C, EquidistantCoordinates<ctype,dim> >,
int> = 0>
912 std::array<
int, std::size_t{dim}> s,
913 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>{0ULL},
917 : ccobj(
comm), _torus(
comm,tag,s,
overlap,partitioner), leafIndexSet_(*
this),
918 _L(L), _periodic(periodic), _coarseSize(s), _overlap(
overlap),
919 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
924 std::fill(o.begin(), o.end(), 0);
925 iTupel o_interior(o);
926 iTupel s_interior(s);
932 for (
int i=0; i<dim; i++)
935 int toosmall = (s_interior[i] < 2*
overlap) &&
936 (periodic[i] || (s_interior[i] != s[i]));
939 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
941 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
942 " Note that this also holds for DOFs on subdomain boundaries."
943 " Increase grid elements or decrease overlap accordingly.");
947 iTupel s_overlap(s_interior);
948 for (
int i=0; i<dim; i++)
950 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
952 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
956 FieldVector<ctype,dim> upperRightWithOverlap;
958 for (
int i=0; i<dim; i++)
959 upperRightWithOverlap[i] = (L[i] / s[i]) * s_overlap[i];
962 EquidistantCoordinates<ctype,dim> cc(upperRightWithOverlap,s_overlap);
979 template<
class C = Coordinates,
980 typename std::enable_if_t< std::is_same_v<C, EquidistantOffsetCoordinates<ctype,dim> >,
int> = 0>
982 Dune::FieldVector<ctype, dim> upperright,
983 std::array<
int, std::size_t{dim}> s,
984 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>(0ULL),
988 : ccobj(
comm), _torus(
comm,tag,s,
overlap,partitioner), leafIndexSet_(*
this),
989 _L(upperright - lowerleft),
990 _periodic(periodic), _coarseSize(s), _overlap(
overlap),
991 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
996 std::fill(o.begin(), o.end(), 0);
997 iTupel o_interior(o);
998 iTupel s_interior(s);
1004 for (
int i=0; i<dim; i++)
1007 int toosmall = (s_interior[i] < 2*
overlap) &&
1008 (periodic[i] || (s_interior[i] != s[i]));
1011 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
1013 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1014 " Note that this also holds for DOFs on subdomain boundaries."
1015 " Increase grid elements or decrease overlap accordingly.");
1019 iTupel s_overlap(s_interior);
1020 for (
int i=0; i<dim; i++)
1022 if ((o_interior[i] -
overlap > 0) || (periodic[i]))
1024 if ((o_interior[i] + s_interior[i] +
overlap <= _coarseSize[i]) || (periodic[i]))
1028 FieldVector<ctype,dim> upperRightWithOverlap;
1029 for (
int i=0; i<dim; i++)
1030 upperRightWithOverlap[i] = lowerleft[i]
1031 + s_overlap[i] * (upperright[i]-lowerleft[i]) / s[i];
1033 EquidistantOffsetCoordinates<ctype,dim> cc(lowerleft,upperRightWithOverlap,s_overlap);
1048 template<
class C = Coordinates,
1049 typename std::enable_if_t< std::is_same_v<C, TensorProductCoordinates<ctype,dim> >,
int> = 0>
1050 YaspGrid (std::array<std::vector<ctype>, std::size_t{dim}> coords,
1051 std::bitset<std::size_t{dim}> periodic = std::bitset<std::size_t{dim}>(0ULL),
1055 : ccobj(
comm), _torus(
comm,tag,Dune::Yasp::sizeArray<dim>(coords),
overlap,partitioner),
1056 leafIndexSet_(*
this), _periodic(periodic), _overlap(
overlap),
1057 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
1060 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1065 for (
int i=0; i<dim; i++) {
1066 _coarseSize[i] = coords[i].size() - 1;
1067 _L[i] = coords[i][_coarseSize[i]] - coords[i][0];
1071 std::fill(o.begin(), o.end(), 0);
1072 iTupel o_interior(o);
1073 iTupel s_interior(_coarseSize);
1075 _torus.
partition(_torus.
rank(),o,_coarseSize,o_interior,s_interior);
1079 for (
int i=0; i<dim; i++)
1082 int toosmall = (s_interior[i] < 2*
overlap) &&
1083 (periodic[i] || (s_interior[i] != _coarseSize[i]));
1086 MPI_Allreduce(&toosmall, &global, 1, MPI_INT, MPI_LOR,
comm);
1088 DUNE_THROW(
Dune::GridError,
"YaspGrid does not support degrees of freedom shared by more than immediately neighboring subdomains."
1089 " Note that this also holds for DOFs on subdomain boundaries."
1090 " Increase grid elements or decrease overlap accordingly.");
1095 std::array<std::vector<ctype>,dim> newcoords;
1096 std::array<int, dim> offset(o_interior);
1099 for (
int i=0; i<dim; ++i)
1102 typename std::vector<ctype>::iterator
begin = coords[i].begin() + o_interior[i];
1103 typename std::vector<ctype>::iterator
end =
begin + s_interior[i] + 1;
1107 if (o_interior[i] -
overlap > 0)
1112 if (o_interior[i] + s_interior[i] +
overlap < _coarseSize[i])
1121 if ((periodic[i]) && (o_interior[i] + s_interior[i] +
overlap >= _coarseSize[i]))
1124 typename std::vector<ctype>::iterator it = coords[i].begin();
1126 newcoords[i].push_back(newcoords[i].back() - *it + *(++it));
1129 if ((periodic[i]) && (o_interior[i] -
overlap <= 0))
1134 typename std::vector<ctype>::iterator it = coords[i].end() - 1;
1136 newcoords[i].insert(newcoords[i].
begin(), newcoords[i].
front() - *it + *(--it));
1140 TensorProductCoordinates<ctype,dim> cc(newcoords, offset);
1163 YaspGrid (std::array<std::vector<ctype>, std::size_t{dim}> coords,
1164 std::bitset<std::size_t{dim}> periodic,
1167 std::array<int,dim> coarseSize,
1169 : ccobj(
comm), _torus(
comm,tag,coarseSize,
overlap,partitioner), leafIndexSet_(*
this),
1170 _periodic(periodic), _coarseSize(coarseSize), _overlap(
overlap),
1171 keep_ovlp(
true), adaptRefCount(0), adaptActive(
false)
1174 static_assert(std::is_same<Coordinates,TensorProductCoordinates<ctype,dim> >::value,
1175 "YaspGrid coordinate container template parameter and given constructor values do not match!");
1178 DUNE_THROW(
Dune::GridError,
"Setup of a tensorproduct grid requires monotonous sequences of coordinates.");
1180 for (
int i=0; i<dim; i++)
1181 _L[i] = coords[i][coords[i].
size() - 1] - coords[i][0];
1185 std::array<int,dim> o;
1186 std::fill(o.begin(), o.end(), 0);
1187 std::array<int,dim> o_interior(o);
1188 std::array<int,dim> s_interior(coarseSize);
1190 _torus.
partition(_torus.
rank(),o,coarseSize,o_interior,s_interior);
1193 std::array<int,dim> offset(o_interior);
1194 for (
int i=0; i<dim; i++)
1195 if ((periodic[i]) || (o_interior[i] > 0))
1198 TensorProductCoordinates<ctype,dim> cc(coords, offset);
1207 friend struct BackupRestoreFacility<
YaspGrid<dim,Coordinates> >;
1219 return _levels.size()-1;
1227 "Coarsening " << -refCount <<
" levels requested!");
1230 for (
int k=refCount; k<0; k++)
1234 _levels.back() = empty;
1238 indexsets.pop_back();
1242 for (
int k=0; k<refCount; k++)
1245 YGridLevel& cg = _levels[
maxLevel()];
1247 std::bitset<dim> ovlp_low(0ULL), ovlp_up(0ULL);
1248 for (
int i=0; i<dim; i++)
1250 if (cg.overlap[0].dataBegin()->origin(i) > 0 || _periodic[i])
1252 if (cg.overlap[0].dataBegin()->max(i) + 1 <
globalSize(i) || _periodic[i])
1256 Coordinates newcont(cg.coords.refine(ovlp_low, ovlp_up, cg.overlapSize, keep_ovlp));
1258 int overlap = (keep_ovlp) ? 2*cg.overlapSize : cg.overlapSize;
1262 for (
int i=0; i<dim; i++)
1263 o_interior[i] = 2*cg.interior[0].dataBegin()->origin(i);
1266 _levels.resize(_levels.size() + 1);
1279 keep_ovlp = keepPhysicalOverlap;
1293 bool mark(
int refCount,
const typename Traits::template Codim<0>::Entity & e )
1295 assert(adaptActive ==
false);
1296 if (e.level() !=
maxLevel())
return false;
1297 adaptRefCount =
std::max(adaptRefCount, refCount);
1307 int getMark (
const typename Traits::template Codim<0>::Entity &e )
const
1309 return ( e.level() ==
maxLevel() ) ? adaptRefCount : 0;
1316 return (adaptRefCount > 0);
1323 adaptRefCount =
comm().max(adaptRefCount);
1324 return (adaptRefCount < 0);
1330 adaptActive =
false;
1335 template<
int cd, PartitionIteratorType pitype>
1336 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator
lbegin (
int level)
const
1338 return levelbegin<cd,pitype>(level);
1342 template<
int cd, PartitionIteratorType pitype>
1343 typename Traits::template Codim<cd>::template Partition<pitype>::LevelIterator
lend (
int level)
const
1345 return levelend<cd,pitype>(level);
1350 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator
lbegin (
int level)
const
1352 return levelbegin<cd,All_Partition>(level);
1357 typename Traits::template Codim<cd>::template Partition<All_Partition>::LevelIterator
lend (
int level)
const
1359 return levelend<cd,All_Partition>(level);
1363 template<
int cd, PartitionIteratorType pitype>
1364 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator
leafbegin ()
const
1366 return levelbegin<cd,pitype>(
maxLevel());
1370 template<
int cd, PartitionIteratorType pitype>
1371 typename Traits::template Codim<cd>::template Partition<pitype>::LeafIterator
leafend ()
const
1373 return levelend<cd,pitype>(
maxLevel());
1378 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator
leafbegin ()
const
1380 return levelbegin<cd,All_Partition>(
maxLevel());
1385 typename Traits::template Codim<cd>::template Partition<All_Partition>::LeafIterator
leafend ()
const
1387 return levelend<cd,All_Partition>(
maxLevel());
1391 template <
typename Seed>
1392 typename Traits::template Codim<Seed::codimension>::Entity
1395 const int codim = Seed::codimension;
1398 typedef typename Traits::template Codim<Seed::codimension>::Entity
Entity;
1402 return Entity(EntityImp(g,YIterator(g->overlapfront[codim],seed.impl().coord(),seed.impl().offset())));
1409 return g->overlapSize;
1416 return g->overlapSize;
1420 int ghostSize ([[maybe_unused]]
int level, [[maybe_unused]]
int codim)
const
1432 int size (
int level,
int codim)
const
1438 typedef typename std::array<YGridComponent<Coordinates>, Dune::power(2,dim)>::iterator DAI;
1439 for (DAI it = g->overlapfront[codim].dataBegin(); it != g->overlapfront[codim].dataEnd(); ++it)
1440 count += it->totalsize();
1454 return (type.isCube()) ?
size(level,dim-type.dim()) : 0;
1478 template<
class DataHandleImp,
class DataType>
1481 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,level);
1488 template<
class DataHandleImp,
class DataType>
1491 YaspCommunicateMeta<dim,dim>::comm(*
this,data,iftype,dir,this->
maxLevel());
1498 template<
class DataHandle,
int codim>
1502 if (!data.contains(dim,codim))
return;
1505 typedef typename DataHandle::DataType DataType;
1516 sendlist = &g->send_interiorborder_interiorborder[codim];
1517 recvlist = &g->recv_interiorborder_interiorborder[codim];
1521 sendlist = &g->send_interiorborder_overlapfront[codim];
1522 recvlist = &g->recv_overlapfront_interiorborder[codim];
1526 sendlist = &g->send_overlap_overlapfront[codim];
1527 recvlist = &g->recv_overlapfront_overlap[codim];
1531 sendlist = &g->send_overlapfront_overlapfront[codim];
1532 recvlist = &g->recv_overlapfront_overlapfront[codim];
1542 std::vector<int> send_size(sendlist->
size(),-1);
1543 std::vector<int> recv_size(recvlist->
size(),-1);
1544 std::vector<size_t*> send_sizes(sendlist->
size(),
static_cast<size_t*
>(0));
1545 std::vector<size_t*> recv_sizes(recvlist->
size(),
static_cast<size_t*
>(0));
1550 if (data.fixedSize(dim,codim))
1554 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1556 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1558 send_size[cnt] = is->grid.totalsize() * data.size(*it);
1562 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1564 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1566 recv_size[cnt] = is->grid.totalsize() * data.size(*it);
1574 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1577 size_t *buf =
new size_t[is->grid.totalsize()];
1578 send_sizes[cnt] = buf;
1581 int i=0;
size_t n=0;
1582 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1584 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1586 for ( ; it!=itend; ++it)
1588 buf[i] = data.size(*it);
1597 torus().
send(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1603 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1606 size_t *buf =
new size_t[is->grid.totalsize()];
1607 recv_sizes[cnt] = buf;
1610 torus().
recv(is->rank,buf,is->grid.totalsize()*
sizeof(
size_t));
1619 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1621 delete[] send_sizes[cnt];
1622 send_sizes[cnt] = 0;
1628 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1631 size_t *buf = recv_sizes[cnt];
1635 for (
int i=0; i<is->grid.totalsize(); ++i)
1646 std::vector<DataType*> sends(sendlist->
size(),
static_cast<DataType*
>(0));
1648 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1651 DataType *buf =
new DataType[send_size[cnt]];
1657 MessageBuffer<DataType> mb(buf);
1660 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1662 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1664 for ( ; it!=itend; ++it)
1665 data.gather(mb,*it);
1668 torus().
send(is->rank,buf,send_size[cnt]*
sizeof(DataType));
1673 std::vector<DataType*> recvs(recvlist->
size(),
static_cast<DataType*
>(0));
1675 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1678 DataType *buf =
new DataType[recv_size[cnt]];
1684 torus().
recv(is->rank,buf,recv_size[cnt]*
sizeof(DataType));
1693 for (ListIt is=sendlist->
begin(); is!=sendlist->
end(); ++is)
1695 delete[] sends[cnt];
1702 for (ListIt is=recvlist->
begin(); is!=recvlist->
end(); ++is)
1705 DataType *buf = recvs[cnt];
1708 MessageBuffer<DataType> mb(buf);
1711 if (data.fixedSize(dim,codim))
1713 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1715 size_t n=data.size(*it);
1716 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1718 for ( ; it!=itend; ++it)
1719 data.scatter(mb,*it,n);
1724 size_t *sbuf = recv_sizes[cnt];
1725 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1727 typename Traits::template Codim<codim>::template Partition<All_Partition>::LevelIterator
1729 for ( ; it!=itend; ++it)
1730 data.scatter(mb,*it,sbuf[i++]);
1743 return theglobalidset;
1748 return theglobalidset;
1753 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1754 return *(indexsets[level]);
1759 return leafIndexSet_;
1782 friend class
Dune::YaspEntity<0, dim, const Dune::YaspGrid<dim, Coordinates> >;
1784 template<int codim_, int dim_, class GridImp_, template<int,int,class> class EntityImp_>
1788 class MessageBuffer {
1791 MessageBuffer (DT *p)
1800 void write (
const Y& data)
1802 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1808 void read (Y& data)
const
1810 static_assert(( std::is_same<DT,Y>::value ),
"DataType mismatch");
1821 template<
int cd, PartitionIteratorType pitype>
1825 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1836 return levelend <cd, pitype> (level);
1838 DUNE_THROW(
GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1842 template<
int cd, PartitionIteratorType pitype>
1843 YaspLevelIterator<cd,pitype,GridImp> levelend (
int level)
const
1846 if (level<0 || level>
maxLevel()) DUNE_THROW(RangeError,
"level out of range");
1849 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interior[cd].end());
1851 return YaspLevelIterator<cd,pitype,GridImp>(g,g->interiorborder[cd].end());
1853 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlap[cd].end());
1855 return YaspLevelIterator<cd,pitype,GridImp>(g,g->overlapfront[cd].end());
1857 DUNE_THROW(GridError,
"YaspLevelIterator with this codim or partition type not implemented");
1862 Torus<Communication,dim> _torus;
1864 std::vector< std::shared_ptr< YaspIndexSet<const YaspGrid<dim,Coordinates>,
false > > > indexsets;
1865 YaspIndexSet<const YaspGrid<dim,Coordinates>,
true> leafIndexSet_;
1866 YaspGlobalIdSet<const YaspGrid<dim,Coordinates> > theglobalidset;
1868 Dune::FieldVector<ctype, dim> _L;
1870 std::bitset<dim> _periodic;
1872 ReservedVector<YGridLevel,32> _levels;
1879#if __cpp_deduction_guides >= 201611
1881 template<
typename ctype,
int dim>
1882 YaspGrid(FieldVector<ctype, dim>,
1883 std::array<
int, std::size_t{dim}>,
1884 std::bitset<std::size_t{dim}> = std::bitset<std::size_t{dim}>{0ULL},
1887 const YLoadBalance<dim>* = YaspGrid< dim, EquidistantCoordinates<ctype, dim> >::defaultLoadbalancer())
1888 -> YaspGrid< dim, EquidistantCoordinates<ctype, dim> >;
1890 template<
typename ctype,
int dim>
1891 YaspGrid(FieldVector<ctype, dim>,
1892 FieldVector<ctype, dim>,
1893 std::array<
int, std::size_t{dim}>,
1894 std::bitset<std::size_t{dim}> = std::bitset<std::size_t{dim}>{0ULL},
1897 const YLoadBalance<dim>* = YaspGrid< dim, EquidistantOffsetCoordinates<ctype, dim> >::defaultLoadbalancer())
1898 -> YaspGrid< dim, EquidistantOffsetCoordinates<ctype, dim> >;
1900 template<
typename ctype, std::
size_t dim>
1901 YaspGrid(std::array<std::vector<ctype>, dim>,
1902 std::bitset<dim> = std::bitset<dim>{0ULL},
1905 const YLoadBalance<
int{dim}>* = YaspGrid<
int{dim}, TensorProductCoordinates<ctype,
int{dim}> >::defaultLoadbalancer())
1906 -> YaspGrid<
int{dim}, TensorProductCoordinates<ctype,
int{dim}> >;
1910 template <
int d,
class CC>
1915 s <<
"[" << rank <<
"]:" <<
" YaspGrid maxlevel=" << grid.
maxLevel() << std::endl;
1917 s <<
"Printing the torus: " <<std::endl;
1918 s << grid.
torus() << std::endl;
1922 s <<
"[" << rank <<
"]: " << std::endl;
1923 s <<
"[" << rank <<
"]: " <<
"==========================================" << std::endl;
1924 s <<
"[" << rank <<
"]: " <<
"level=" << g->level() << std::endl;
1926 for (
int codim = 0; codim < d + 1; ++codim)
1928 s <<
"[" << rank <<
"]: " <<
"overlapfront[" << codim <<
"]: " << g->overlapfront[codim] << std::endl;
1929 s <<
"[" << rank <<
"]: " <<
"overlap[" << codim <<
"]: " << g->overlap[codim] << std::endl;
1930 s <<
"[" << rank <<
"]: " <<
"interiorborder[" << codim <<
"]: " << g->interiorborder[codim] << std::endl;
1931 s <<
"[" << rank <<
"]: " <<
"interior[" << codim <<
"]: " << g->interior[codim] << std::endl;
1934 for (I i=g->send_overlapfront_overlapfront[codim].begin();
1935 i!=g->send_overlapfront_overlapfront[codim].end(); ++i)
1936 s <<
"[" << rank <<
"]: " <<
" s_of_of[" << codim <<
"] to rank "
1937 << i->rank <<
" " << i->grid << std::endl;
1939 for (I i=g->recv_overlapfront_overlapfront[codim].begin();
1940 i!=g->recv_overlapfront_overlapfront[codim].end(); ++i)
1941 s <<
"[" << rank <<
"]: " <<
" r_of_of[" << codim <<
"] to rank "
1942 << i->rank <<
" " << i->grid << std::endl;
1944 for (I i=g->send_overlap_overlapfront[codim].begin();
1945 i!=g->send_overlap_overlapfront[codim].end(); ++i)
1946 s <<
"[" << rank <<
"]: " <<
" s_o_of[" << codim <<
"] to rank "
1947 << i->rank <<
" " << i->grid << std::endl;
1949 for (I i=g->recv_overlapfront_overlap[codim].begin();
1950 i!=g->recv_overlapfront_overlap[codim].end(); ++i)
1951 s <<
"[" << rank <<
"]: " <<
" r_of_o[" << codim <<
"] to rank "
1952 << i->rank <<
" " << i->grid << std::endl;
1954 for (I i=g->send_interiorborder_interiorborder[codim].begin();
1955 i!=g->send_interiorborder_interiorborder[codim].end(); ++i)
1956 s <<
"[" << rank <<
"]: " <<
" s_ib_ib[" << codim <<
"] to rank "
1957 << i->rank <<
" " << i->grid << std::endl;
1959 for (I i=g->recv_interiorborder_interiorborder[codim].begin();
1960 i!=g->recv_interiorborder_interiorborder[codim].end(); ++i)
1961 s <<
"[" << rank <<
"]: " <<
" r_ib_ib[" << codim <<
"] to rank "
1962 << i->rank <<
" " << i->grid << std::endl;
1964 for (I i=g->send_interiorborder_overlapfront[codim].begin();
1965 i!=g->send_interiorborder_overlapfront[codim].end(); ++i)
1966 s <<
"[" << rank <<
"]: " <<
" s_ib_of[" << codim <<
"] to rank "
1967 << i->rank <<
" " << i->grid << std::endl;
1969 for (I i=g->recv_overlapfront_interiorborder[codim].begin();
1970 i!=g->recv_overlapfront_interiorborder[codim].end(); ++i)
1971 s <<
"[" << rank <<
"]: " <<
" r_of_ib[" << codim <<
"] to rank "
1972 << i->rank <<
" " << i->grid << std::endl;
1981 namespace Capabilities
1991 template<
int dim,
class Coordinates>
1994 static const bool v =
true;
2000 template<
int dim,
class Coordinates>
2003 static const bool v =
true;
2004 static const unsigned int topologyId = GeometryTypes::cube(dim).id();
2010 template<
int dim,
class Coordinates>
2013 static const bool v =
true;
2019 template<
int dim,
class Coordinates,
int codim>
2022 static const bool v =
true;
2029 template<
int dim,
class Coordinates,
int codim>
2032 static const bool v =
true;
2038 template<
int dim,
int codim,
class Coordinates>
2041 static const bool v =
true;
2047 template<
int dim,
class Coordinates>
2050 static const bool v =
true;
2056 template<
int dim,
class Coordinates>
2059 static const bool v =
true;
2065 template<
int dim,
class Coordinates>
2068 static const bool v =
true;
This provides a YGrid, the elemental component of the yaspgrid implementation.
Specialization of the PersistentContainer for YaspGrid.
The YaspLevelIterator class.
The YaspIntersectionIterator class.
The YaspIntersection class.
level-wise, non-persistent, consecutive indices for YaspGrid
The YaspGeometry class and its specializations.
The YaspEntitySeed class.
the YaspEntity class and its specializations
This file provides the infrastructure for toroidal communication in YaspGrid.
Specialization of the StructuredGridFactory class for YaspGrid.
This provides container classes for the coordinates to be used in YaspGrid Upon implementation of the...
Provides base classes for index and id sets.
Describes the parallel communication interface class for MessageBuffers and DataHandles.
std::ostream & operator<<(std::ostream &out, const PartitionType &type)
write a PartitionType to a stream
Definition: gridenums.hh:72
CommunicationDirection
Define a type for communication direction parameter.
Definition: gridenums.hh:170
InterfaceType
Parameter to be used for the communication functions.
Definition: gridenums.hh:86
@ All_Partition
all entities
Definition: gridenums.hh:141
@ Interior_Partition
only interior entities
Definition: gridenums.hh:137
@ InteriorBorder_Partition
interior and border entities
Definition: gridenums.hh:138
@ Overlap_Partition
interior, border, and overlap entities
Definition: gridenums.hh:139
@ Ghost_Partition
only ghost entities
Definition: gridenums.hh:142
@ BackwardCommunication
reverse communication direction
Definition: gridenums.hh:172
@ InteriorBorder_All_Interface
send interior and border, receive all entities
Definition: gridenums.hh:88
@ All_All_Interface
send all and receive all entities
Definition: gridenums.hh:91
@ Overlap_All_Interface
send overlap, receive all entities
Definition: gridenums.hh:90
@ Overlap_OverlapFront_Interface
send overlap, receive overlap and front entities
Definition: gridenums.hh:89
@ InteriorBorder_InteriorBorder_Interface
send/receive interior and border entities
Definition: gridenums.hh:87
void swap(Dune::PersistentContainer< G, T > &a, Dune::PersistentContainer< G, T > &b)
Definition: utility/persistentcontainer.hh:83
Include standard header files.
Definition: agrid.hh:60
const int yaspgrid_level_bits
Definition: yaspgrid.hh:48
Communication< MPI_Comm > YaspCommunication
Definition: yaspgrid.hh:85
const int yaspgrid_dim_bits
Definition: yaspgrid.hh:47
int max(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:337
constexpr Overlap overlap
PartitionSet for the overlap partition.
Definition: partitionset.hh:278
constexpr Front front
PartitionSet for the front partition.
Definition: partitionset.hh:281
constexpr Interior interior
PartitionSet for the interior partition.
Definition: partitionset.hh:272
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:132
bool checkIfMonotonous(const std::array< std::vector< ctype >, dim > &coords)
Definition: coordinates.hh:372
Intersection of a mesh entity of codimension 0 ("element") with a "neighboring" element or with the d...
Definition: common/intersection.hh:164
Specialize with 'true' for if the codimension 0 entity of the grid has only one possible geometry typ...
Definition: common/capabilities.hh:27
static const bool v
Definition: common/capabilities.hh:28
static const unsigned int topologyId
Definition: common/capabilities.hh:31
Specialize with 'true' if the grid is a Cartesian grid. Cartesian grids satisfy the following propert...
Definition: common/capabilities.hh:48
static const bool v
Definition: common/capabilities.hh:50
Specialize with 'true' for all codims that a grid implements entities for. (default=false)
Definition: common/capabilities.hh:58
static const bool v
Definition: common/capabilities.hh:59
specialize with 'true' for all codims that a grid provides an iterator for (default=hasEntity<codim>:...
Definition: common/capabilities.hh:74
static const bool v
Definition: common/capabilities.hh:75
specialize with 'true' for all codims that a grid can communicate data on (default=false)
Definition: common/capabilities.hh:97
static const bool v
Definition: common/capabilities.hh:98
Specialize with 'true' if implementation guarantees conforming level grids. (default=false)
Definition: common/capabilities.hh:106
static const bool v
Definition: common/capabilities.hh:107
Specialize with 'true' if implementation guarantees a conforming leaf grid. (default=false)
Definition: common/capabilities.hh:115
static const bool v
Definition: common/capabilities.hh:116
Specialize with 'true' if implementation provides backup and restore facilities. (default=false)
Definition: common/capabilities.hh:124
static const bool v
Definition: common/capabilities.hh:125
Specialize with 'true' if the grid implementation is thread safe, while it is not modified....
Definition: common/capabilities.hh:169
static const bool v
Definition: common/capabilities.hh:170
CommDataHandleIF describes the features of a data handle for communication in parallel runs using the...
Definition: datahandleif.hh:78
Definition: defaultgridview.hh:26
Definition: defaultgridview.hh:219
Wrapper class for entities.
Definition: common/entity.hh:66
Base class for exceptions in Dune grid modules.
Definition: exceptions.hh:20
Definition: common/grid.hh:862
Index Set Interface base class.
Definition: indexidset.hh:78
Id Set Interface.
Definition: indexidset.hh:452
detected_or_fallback_t< DeprecatedCollectiveCommunication_t, Communication_t, typename GridFamily::Traits > Communication
A type that is a model of Dune::Communication. It provides a portable way for communication on the se...
Definition: common/grid.hh:525
A traits struct that collects all associated types of one grid model.
Definition: common/grid.hh:995
[ provides Dune::Grid ]
Definition: yaspgrid.hh:163
YaspGridFamily< dim, Coordinates > GridFamily
the GridFamily of this grid
Definition: yaspgrid.hh:722
YaspIndexSet< YaspGrid< dim, Coordinates >, false > LevelIndexSetType
Definition: yaspgrid.hh:727
int ghostSize(int level, int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1420
const Torus< Communication, dim > & torus() const
return reference to torus
Definition: yaspgrid.hh:244
typename Base::Communication Communication
Definition: yaspgrid.hh:178
const Traits::LeafIndexSet & leafIndexSet() const
Definition: yaspgrid.hh:1757
void init()
Definition: yaspgrid.hh:684
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lbegin(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1350
int size(int level, GeometryType type) const
number of entities per level and geometry type in this process
Definition: yaspgrid.hh:1452
Traits::template Codim< Seed::codimension >::Entity entity(const Seed &seed) const
Definition: yaspgrid.hh:1393
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1378
const Traits::GlobalIdSet & globalIdSet() const
Definition: yaspgrid.hh:1741
const Traits::LocalIdSet & localIdSet() const
Definition: yaspgrid.hh:1746
YaspGridFamily< dim, Coordinates >::Traits Traits
Definition: yaspgrid.hh:724
bool getRefineOption() const
Definition: yaspgrid.hh:282
void globalRefine(int refCount)
refine the grid refCount times.
Definition: yaspgrid.hh:1223
int getMark(const typename Traits::template Codim< 0 >::Entity &e) const
returns adaptation mark for given entity
Definition: yaspgrid.hh:1307
void boundarysegmentssize()
Definition: yaspgrid.hh:690
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, Communication comm=Communication(), const Yasp::Partitioning< dim > *partitioner=defaultPartitioner())
Standard constructor for a tensorproduct YaspGrid.
Definition: yaspgrid.hh:1050
YaspIndexSet< YaspGrid< dim, Coordinates >, true > LeafIndexSetType
Definition: yaspgrid.hh:728
int ghostSize(int codim) const
return size (= distance in graph) of ghost region
Definition: yaspgrid.hh:1426
typename Base::Communication CollectiveCommunicationType
Definition: yaspgrid.hh:179
int globalSize(int i) const
return number of cells on finest level in given direction on all processors
Definition: yaspgrid.hh:250
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, Communication comm=Communication(), const Yasp::Partitioning< dim > *partitioner=defaultPartitioner())
Definition: yaspgrid.hh:981
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1371
void postAdapt()
clean up some markers
Definition: yaspgrid.hh:1328
const Dune::FieldVector< ctype, dim > & domainSize() const
returns the size of the physical domain
Definition: yaspgrid.hh:1470
YGridLevelIterator end() const
return iterator pointing to one past the finest level
Definition: yaspgrid.hh:305
friend class Entity
Definition: yaspgrid.hh:1785
int size(GeometryType type) const
number of leaf entities per geometry type in this process
Definition: yaspgrid.hh:1458
int maxLevel() const
Definition: yaspgrid.hh:1217
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, Communication comm=Communication(), const Yasp::Partitioning< dim > *partitioner=defaultPartitioner())
Definition: yaspgrid.hh:911
Traits::template Codim< cd >::template Partition< All_Partition >::LeafIterator leafend() const
return LeafIterator which points behind the last entity in maxLevel
Definition: yaspgrid.hh:1385
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1479
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:568
Traits::template Codim< cd >::template Partition< pitype >::LeafIterator leafbegin() const
return LeafIterator which points to the first entity in maxLevel
Definition: yaspgrid.hh:1364
bool preAdapt()
returns true, if the grid will be coarsened
Definition: yaspgrid.hh:1320
iTupel levelSize(int l) const
return size vector of the grid (in cells) on level l
Definition: yaspgrid.hh:268
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:1293
int overlapSize(int level, int codim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1406
void communicate(CommDataHandleIF< DataHandleImp, DataType > &data, InterfaceType iftype, CommunicationDirection dir) const
Definition: yaspgrid.hh:1489
int size(int codim) const
number of leaf entities per codim in this process
Definition: yaspgrid.hh:1446
bool isPeriodic(int i) const
return whether the grid is periodic in direction i
Definition: yaspgrid.hh:277
const Traits::LevelIndexSet & levelIndexSet(int level) const
Definition: yaspgrid.hh:1751
void refineOptions(bool keepPhysicalOverlap)
set options for refinement
Definition: yaspgrid.hh:1277
int size(int level, int codim) const
number of entities per level and codim in this process
Definition: yaspgrid.hh:1432
YaspGlobalIdSet< YaspGrid< dim, Coordinates > > GlobalIdSetType
Definition: yaspgrid.hh:729
YaspGrid(const Coordinates &coordinates, std::bitset< dim > periodic=std::bitset< dim >(0ULL), int overlap=1, Communication comm=Communication(), const Yasp::Partitioning< dim > *partitioner=defaultPartitioner())
Definition: yaspgrid.hh:738
bool adapt()
map adapt to global refine
Definition: yaspgrid.hh:1313
Traits::template Codim< cd >::template Partition< All_Partition >::LevelIterator lend(int level) const
version without second template parameter for convenience
Definition: yaspgrid.hh:1357
ReservedVector< YGridLevel, 32 >::const_iterator YGridLevelIterator
Iterator over the grid levels.
Definition: yaspgrid.hh:288
size_t numBoundarySegments() const
returns the number of boundary segments within the macro grid
Definition: yaspgrid.hh:1464
void communicateCodim(DataHandle &data, InterfaceType iftype, CommunicationDirection dir, int level) const
Definition: yaspgrid.hh:1499
void makelevel(const Coordinates &coords, std::bitset< dim > periodic, iTupel o_interior, int overlap)
Make a new YGridLevel structure.
Definition: yaspgrid.hh:333
bigunsignedint< dim *yaspgrid_dim_bits+yaspgrid_level_bits+dim > PersistentIndexType
Definition: yaspgrid.hh:719
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:1343
static const YLoadBalanceDefault< dim > * defaultLoadbalancer()
Definition: yaspgrid.hh:312
Traits::template Codim< cd >::template Partition< pitype >::LevelIterator lbegin(int level) const
one past the end on this level
Definition: yaspgrid.hh:1336
static const Yasp::Partitioning< dim > * defaultPartitioner()
Definition: yaspgrid.hh:319
iTupel globalSize() const
return number of cells on finest level on all processors
Definition: yaspgrid.hh:256
const Communication & comm() const
return a communication object
Definition: yaspgrid.hh:1764
int overlapSize(int odim) const
return size (= distance in graph) of overlap region
Definition: yaspgrid.hh:1413
const YaspGrid< dim, Coordinates > GridImp
Definition: yaspgrid.hh:682
YGridLevelIterator begin(int i) const
return iterator pointing to given level
Definition: yaspgrid.hh:297
Coordinates::ctype ctype
Type used for coordinates.
Definition: yaspgrid.hh:176
int levelSize(int l, int i) const
return size of the grid (in cells) on level l in direction i
Definition: yaspgrid.hh:262
YGridLevelIterator begin() const
return iterator pointing to coarsest level
Definition: yaspgrid.hh:291
The general version that handles all codimensions but 0 and dim.
Definition: yaspgridgeometry.hh:31
Definition: yaspgridentity.hh:268
Describes the minimal information necessary to create a fully functional YaspEntity.
Definition: yaspgridentityseed.hh:18
Iterates over entities of one grid level.
Definition: yaspgridleveliterator.hh:19
YaspIntersectionIterator enables iteration over intersections with neighboring codim 0 entities.
Definition: yaspgridintersectioniterator.hh:22
YaspIntersection provides data about intersection with neighboring codim 0 entities.
Definition: yaspgridintersection.hh:22
YaspHierarchicIterator enables iteration over son entities of codim 0.
Definition: yaspgridhierarchiciterator.hh:20
Implementation of Level- and LeafIndexSets for YaspGrid.
Definition: yaspgridindexsets.hh:25
persistent, globally unique Ids
Definition: yaspgrididset.hh:25
Definition: yaspgridpersistentcontainer.hh:35
Definition: yaspgrid.hh:92
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:115
YaspCommunication CCType
Definition: yaspgrid.hh:93
Container for equidistant coordinates in a YaspGrid.
Definition: coordinates.hh:29
Container for equidistant coordinates in a YaspGrid with non-trivial origin.
Definition: coordinates.hh:131
Coordinate container for a tensor product YaspGrid.
Definition: coordinates.hh:245
a base class for the yaspgrid partitioning strategy
Definition: partitioning.hh:39
Definition: partitioning.hh:48
Implement the default load balance strategy of yaspgrid.
Definition: partitioning.hh:206
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:239
iTupel coord() const
return own coordinates
Definition: torus.hh:100
int rank() const
return own rank
Definition: torus.hh:94
const iTupel & dims() const
return dimensions of torus
Definition: torus.hh:112
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:374
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:361
int neighbors() const
return the number of neighbors, which is
Definition: torus.hh:203
ProcListIterator recvend() const
last process in receive list
Definition: torus.hh:355
ProcListIterator sendend() const
end of send list
Definition: torus.hh:343
ProcListIterator sendbegin() const
first process in send list
Definition: torus.hh:337
void exchange() const
exchange messages stored in request buffers; clear request buffers afterwards
Definition: torus.hh:387
ProcListIterator recvbegin() const
first process in receive list
Definition: torus.hh:349
YGridComponent< Coordinates > move(iTupel v) const
return grid moved by the vector v
Definition: ygrid.hh:263
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:271
implements a collection of YGridComponents which form a codimension Entities of given codimension c n...
Definition: ygrid.hh:551
Iterator over a collection o YGrids A YGrid::Iterator is the heart of an entity in YaspGrid.
Definition: ygrid.hh:594
implements a collection of multiple std::deque<Intersection> Intersections with neighboring processor...
Definition: ygrid.hh:823
int size() const
return the size of the container, this is the sum of the sizes of all deques
Definition: ygrid.hh:953
Iterator end() const
return iterator pointing to the end of the container
Definition: ygrid.hh:929
Iterator begin() const
return iterator pointing to the begin of the container
Definition: ygrid.hh:923
type describing an intersection with a neighboring processor
Definition: ygrid.hh:829
A set of traits classes to store static information about grid implementation.
Different resources needed by all grid implementations.