3 #ifndef DUNE_GRID_YASPGRIDINTERSECTION_HH
4 #define DUNE_GRID_YASPGRIDINTERSECTION_HH
18 template<
class Gr
idImp>
19 class YaspIntersection
21 enum { dim=GridImp::dimension };
22 enum { dimworld=GridImp::dimensionworld };
23 typedef typename GridImp::ctype ctype;
25 typedef typename GridImp::Traits::template Codim< 1 >::GeometryImpl GeometryImpl;
26 typedef typename GridImp::Traits::template Codim< 1 >::LocalGeometryImpl LocalGeometryImpl;
32 typedef typename GridImp::YGridLevelIterator
YGLI;
33 typedef typename GridImp::YGrid::Iterator
I;
34 typedef typename GridImp::template Codim<0>::Entity
Entity;
35 typedef typename GridImp::template Codim<1>::Geometry
Geometry;
36 typedef typename GridImp::template Codim<1>::LocalGeometry
LocalGeometry;
41 std::array<int,dim> dist{{0}};
44 dist[_dir] = 1 - 2*_face;
51 dist[_dir] += -1 + 2*_face;
64 if (_inside.
gridlevel()->mg->isPeriodic(_dir))
77 return coord > _inside.
gridlevel()->overlap[0].dataBegin()->min(_dir)
79 coord <= _inside.
gridlevel()->overlap[0].dataBegin()->max(_dir);
106 DUNE_THROW(
GridError,
"called boundarySegmentIndex while boundary() == false");
108 const std::array<int, dim> & size = _inside.
gridlevel()->mg->begin()->overlap[0].dataBegin()->size();
109 const std::array<int, dim> & origin = _inside.
gridlevel()->mg->begin()->overlap[0].dataBegin()->origin();
110 std::array<int, dim> sides;
112 for (
int i=0; i<dim; i++)
115 ((_inside.
gridlevel()->mg->begin()->overlap[0].dataBegin()->origin(i)
117 (_inside.
gridlevel()->mg->begin()->overlap[0].dataBegin()->origin(i) +
118 _inside.
gridlevel()->mg->begin()->overlap[0].dataBegin()->size(i)
120 _inside.
gridlevel()->mg->levelSize(0,i)));
126 for(
int i=0; i<dim; i++)
128 pos[i] = pos[i] / (1<<_inside.
level());
129 pos[i] = pos[i] - origin[i];
132 std::array<int, dim> fsize;
135 for (
int k=0; k<dim; k++)
137 for (
int k=0; k<dim; k++)
138 fsize[k] = vol/size[k];
144 for (
int k=dim-1; k>=0; k--)
146 if (k == _dir)
continue;
147 index += (pos[k]) * localoffset;
148 localoffset *= size[k];
153 for (
int k=0; k<_dir; k++)
154 index += sides[k] * fsize[k];
156 index += _face * (sides[_dir]>1) * fsize[_dir];
163 FieldVector<ctype, dimworld>
outerNormal (
const FieldVector<ctype, dim-1>& local)
const
165 return _faceInfo[_count].normal;
169 FieldVector<ctype, dimworld>
unitOuterNormal (
const FieldVector<ctype, dim-1>& local)
const
171 return _faceInfo[_count].normal;
177 return _faceInfo[_count].normal;
185 FieldVector<ctype, dimworld> n = _faceInfo[_count].normal;
211 std::bitset<dim> shift;
215 Dune::FieldVector<ctype,dimworld> ll, ur;
216 for (
int i=0; i<dimworld; i++)
220 if ((i == _dir) and (_face))
229 if (_inside.
gridlevel()->mg->isPeriodic(i)) {
231 if (coordPeriodic < 0) {
232 auto size = _inside.
gridlevel()->mg->domainSize()[i];
235 }
else if (coordPeriodic + 1 > _inside.
gridlevel()->mg->levelSize(_inside.
gridlevel()->level(),i)) {
236 auto size = _inside.
gridlevel()->mg->domainSize()[i];
243 GeometryImpl _is_global(ll,ur,shift);
250 return GeometryTypes::cube(dim-1);
274 _inside(myself.gridlevel(),
275 myself.transformingsubiterator()),
276 _outside(myself.gridlevel(),
277 myself.transformingsubiterator()),
306 return _count == other._count && _inside.
equals(other._inside);
321 FieldVector<ctype, dimworld> normal;
322 LocalGeometryImpl geom_inside;
323 LocalGeometryImpl geom_outside;
327 static const std::array<faceInfo, 2*GridImp::dimension> _faceInfo;
329 static std::array<faceInfo, 2*dim> initFaceInfo()
331 std::array<faceInfo, 2*dim>
I;
336 I[2*i+1].normal = 0.0;
337 I[2*i].normal[i] = -1.0;
338 I[2*i+1].normal[i] = +1.0;
346 Dune::FieldVector<ctype, dim> ll(0.0);
347 Dune::FieldVector<ctype, dim> ur(1.0);
350 I[2*i].geom_inside = LocalGeometryImpl(ll,ur,s);
351 I[2*i+1].geom_outside = LocalGeometryImpl(ll,ur,s);
356 I[2*i].geom_outside = LocalGeometryImpl(ll,ur,s);
357 I[2*i+1].geom_inside = LocalGeometryImpl(ll,ur,s);
364 template<
class Gr
idImp>
365 const std::array<typename YaspIntersection<GridImp>::faceInfo, 2*GridImp::dimension>
366 YaspIntersection<GridImp>::_faceInfo =
367 YaspIntersection<GridImp>::initFaceInfo();
371 #endif // DUNE_GRID_YASPGRIDINTERSECTION_HH