3.2-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
3dinteractionvolumecontaineradaptive.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_FVMPFAL3D_INTERACTIONVOLUMECONTAINER_ADAPTIVE_HH
25#define DUMUX_FVMPFAL3D_INTERACTIONVOLUMECONTAINER_ADAPTIVE_HH
26
27// dumux environment
32
33namespace Dumux {
34
45template<class TypeTag>
47{
48 friend class FvMpfaL3dInteractionVolumeContainer<TypeTag>;
51
52 enum
53 {
54 dim = GridView::dimension, dimWorld = GridView::dimensionworld
55 };
56
59
60 using ReferenceElements = Dune::ReferenceElements<Scalar, dim>;
61
63
64 using Element = typename GridView::Traits::template Codim<0>::Entity;
65 using Vertex = typename GridView::Traits::template Codim<dim>::Entity;
66 using ElementGeometry = typename Element::Geometry;
67
68 using GlobalPosition = typename ElementGeometry::GlobalCoordinate;
69 using DimMatrix = Dune::FieldMatrix<Scalar, dim, dim>;
70
71 using DimVector = Dune::FieldVector<Scalar, dim>;
72 using FaceVerticesVector = std::vector<Dune::FieldVector<std::set<int>, 2*dim> >;
73
74 enum
75 {
76 pressureEqIdx = Indices::pressureEqIdx,
77 };
78
79 using IndexTranslator = IndexTranslatorAdaptive;
80
81public:
85
86private:
87
88 void storeHangingNodeInteractionVolume(InteractionVolume& interactionVolume, const Vertex& vertex);
89 void storeInnerInteractionVolume(InteractionVolume& interactionVolume, const Vertex& vertex);
90
91protected:
93
94public:
95
104 std::set<int>& faceVerticeIndices(int eIdxGlobal, int fIdx)
105 {
106 return faceVertices_[eIdxGlobal][fIdx];
107 }
108
114 ParentType(problem), problem_(problem)
115 {
116 if (dim != 3)
117 {
118 DUNE_THROW(Dune::NotImplemented, "Dimension not supported!");
119 }
120 }
121
122private:
124
126 Implementation &asImp_()
127 { return *static_cast<Implementation *>(this);}
128
130 const Implementation &asImp_() const
131 { return *static_cast<const Implementation *>(this);}
132
133 Problem& problem_;
134 FaceVerticesVector faceVertices_;
135};
136
160template<class TypeTag>
161void FvMpfaL3dInteractionVolumeContainerAdaptive<TypeTag>::storeInnerInteractionVolume(InteractionVolume& interactionVolume,
162 const Vertex& vertex)
163{
164 // if cells of different grid level appear, the geometric information is constructed going from the
165 // coarsest level to the finest level. This ensures that coarser information ins always overwritten by
166 // finer information.
167 if (!interactionVolume.sameLevel())
168 {
169 // sort the local element indices according to the grid level
170 std::vector < std::vector<int> > levelIdx(8, std::vector<int>(2));
171 for (int i = 0; i < 8; i++)
172 {
173 levelIdx[i][0] = i;
174 levelIdx[i][1] = interactionVolume.getSubVolumeElement(i).level();
175 }
176
177 std::sort(levelIdx.begin(), levelIdx.end(), [](const auto& a, const auto& b) { return (a[1]<b[1]); });
178
179 // Generate and store the geometric information going from the coarsest to the finest level.
180 // For the calculation we take advantage from the fact that the ordering inside the interaction volume
181 // with respect to the DUNE reference element is known due to the storage process of the elements in
182 // FvMpfaL3dInteractionVolumeContainer<TypeTag>::storeSubVolumeElements
183 // (const Element& element, std::vector < std::vector<int> >& elemVertMap)
184 for (int i = 0; i < 8; i++)
185 {
186 int idx = levelIdx[i][0];
187 auto element = interactionVolume.getSubVolumeElement(idx);
188
189 const ElementGeometry& geometry = element.geometry();
190
191 const auto referenceElement = ReferenceElements::general(geometry.type());
192
193 switch (idx)
194 {
195 case 0:
196 {
197 DimVector edgeCoord(geometry.global(referenceElement.position(9, dim - 1)));
198 interactionVolume.setEdgePosition(edgeCoord, 2);
199 edgeCoord = geometry.global(referenceElement.position(3, dim - 1));
200 interactionVolume.setEdgePosition(edgeCoord, 0);
201 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
202 interactionVolume.setEdgePosition(edgeCoord, 5);
203
204 break;
205 }
206 case 1:
207 {
208 DimVector edgeCoord(geometry.global(referenceElement.position(2, dim - 1)));
209 interactionVolume.setEdgePosition(edgeCoord, 0);
210 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
211 interactionVolume.setEdgePosition(edgeCoord, 2);
212 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
213 interactionVolume.setEdgePosition(edgeCoord, 3);
214
215 break;
216 }
217 case 2:
218 {
219 DimVector edgeCoord(geometry.global(referenceElement.position(1, dim - 1)));
220 interactionVolume.setEdgePosition(edgeCoord, 0);
221 edgeCoord = geometry.global(referenceElement.position(9, dim - 1));
222 interactionVolume.setEdgePosition(edgeCoord, 4);
223 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
224 interactionVolume.setEdgePosition(edgeCoord, 5);
225
226 break;
227 }
228 case 3:
229 {
230 DimVector edgeCoord(geometry.global(referenceElement.position(0, dim - 1)));
231 interactionVolume.setEdgePosition(edgeCoord, 0);
232 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
233 interactionVolume.setEdgePosition(edgeCoord, 4);
234 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
235 interactionVolume.setEdgePosition(edgeCoord, 3);
236
237 break;
238 }
239 case 4:
240 {
241 DimVector edgeCoord(geometry.global(referenceElement.position(3, dim - 1)));
242 interactionVolume.setEdgePosition(edgeCoord, 1);
243 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
244 interactionVolume.setEdgePosition(edgeCoord, 2);
245 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
246 interactionVolume.setEdgePosition(edgeCoord, 5);
247
248 break;
249 }
250 case 5:
251 {
252 DimVector edgeCoord(geometry.global(referenceElement.position(2, dim - 1)));
253 interactionVolume.setEdgePosition(edgeCoord, 1);
254 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
255 interactionVolume.setEdgePosition(edgeCoord, 2);
256 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
257 interactionVolume.setEdgePosition(edgeCoord, 3);
258
259 break;
260 }
261 case 6:
262 {
263 DimVector edgeCoord(geometry.global(referenceElement.position(1, dim - 1)));
264 interactionVolume.setEdgePosition(edgeCoord, 1);
265 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
266 interactionVolume.setEdgePosition(edgeCoord, 4);
267 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
268 interactionVolume.setEdgePosition(edgeCoord, 5);
269
270 break;
271 }
272 case 7:
273 {
274 DimVector edgeCoord(geometry.global(referenceElement.position(0, dim - 1)));
275 interactionVolume.setEdgePosition(edgeCoord, 1);
276 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
277 interactionVolume.setEdgePosition(edgeCoord, 4);
278 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
279 interactionVolume.setEdgePosition(edgeCoord, 3);
280
281 break;
282 }
283 }
284 }
285 ParentType::storeInnerInteractionVolume(interactionVolume, vertex, false);
286 }
287 else
288 {
289 ParentType::storeInnerInteractionVolume(interactionVolume, vertex, true);
290 }
291}
292
315template<class TypeTag>
316void FvMpfaL3dInteractionVolumeContainerAdaptive<TypeTag>::storeHangingNodeInteractionVolume(InteractionVolume& interactionVolume,
317 const Vertex& vertex)
318{
319 const DimVector& centerPos = vertex.geometry().center();
320
321 interactionVolume.setCenterPosition(centerPos);
322
323 // sort the local element indices according to the grid level
324 std::vector < std::vector<int> > levelIdx(8, std::vector<int>(2));
325 for (int i = 0; i < 8; i++)
326 {
327 levelIdx[i][0] = i;
328 if (interactionVolume.hasSubVolumeElement(i))
329 levelIdx[i][1] = interactionVolume.getSubVolumeElement(i).level();
330 else
331 levelIdx[i][1] = -1;
332 }
333
334 std::sort(levelIdx.begin(), levelIdx.end(), [](const auto& a, const auto& b) { return (a[1]<b[1]); });
335
336 // Generate and store the geometric information going from the coarsest to the finest level.
337 // For the calculation we take advantage from the fact that the ordering inside the interaction volume
338 // with respect to the DUNE reference element is known due to the storage process of the elements in
339 // FvMpfaL3dInteractionVolumeContainer<TypeTag>::
340 // storeSubVolumeElements(const Element& element, std::vector < std::vector<int> >& elemVertMap)
341 for (int i = 0; i < 8; i++)
342 {
343 if (levelIdx[i][1] < 0)
344 continue;
345
346 int idx = levelIdx[i][0];
347
348 auto element = interactionVolume.getSubVolumeElement(idx);
349
350 const ElementGeometry& geometry = element.geometry();
351
352 const auto referenceElement = ReferenceElements::general(geometry.type());
353
354 switch (idx)
355 {
356 case 0:
357 {
358 DimVector edgeCoord(geometry.global(referenceElement.position(9, dim - 1)));
359 interactionVolume.setEdgePosition(edgeCoord, 2);
360 edgeCoord = geometry.global(referenceElement.position(3, dim - 1));
361 interactionVolume.setEdgePosition(edgeCoord, 0);
362 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
363 interactionVolume.setEdgePosition(edgeCoord, 5);
364
365 break;
366 }
367 case 1:
368 {
369 DimVector edgeCoord(geometry.global(referenceElement.position(2, dim - 1)));
370 interactionVolume.setEdgePosition(edgeCoord, 0);
371 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
372 interactionVolume.setEdgePosition(edgeCoord, 2);
373 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
374 interactionVolume.setEdgePosition(edgeCoord, 3);
375
376 break;
377 }
378 case 2:
379 {
380 DimVector edgeCoord(geometry.global(referenceElement.position(1, dim - 1)));
381 interactionVolume.setEdgePosition(edgeCoord, 0);
382 edgeCoord = geometry.global(referenceElement.position(9, dim - 1));
383 interactionVolume.setEdgePosition(edgeCoord, 4);
384 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
385 interactionVolume.setEdgePosition(edgeCoord, 5);
386
387 break;
388 }
389 case 3:
390 {
391 DimVector edgeCoord(geometry.global(referenceElement.position(0, dim - 1)));
392 interactionVolume.setEdgePosition(edgeCoord, 0);
393 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
394 interactionVolume.setEdgePosition(edgeCoord, 4);
395 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
396 interactionVolume.setEdgePosition(edgeCoord, 3);
397
398 break;
399 }
400 case 4:
401 {
402 DimVector edgeCoord(geometry.global(referenceElement.position(3, dim - 1)));
403 interactionVolume.setEdgePosition(edgeCoord, 1);
404 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
405 interactionVolume.setEdgePosition(edgeCoord, 2);
406 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
407 interactionVolume.setEdgePosition(edgeCoord, 5);
408
409 break;
410 }
411 case 5:
412 {
413 DimVector edgeCoord(geometry.global(referenceElement.position(2, dim - 1)));
414 interactionVolume.setEdgePosition(edgeCoord, 1);
415 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
416 interactionVolume.setEdgePosition(edgeCoord, 2);
417 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
418 interactionVolume.setEdgePosition(edgeCoord, 3);
419
420 break;
421 }
422 case 6:
423 {
424 DimVector edgeCoord(geometry.global(referenceElement.position(1, dim - 1)));
425 interactionVolume.setEdgePosition(edgeCoord, 1);
426 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
427 interactionVolume.setEdgePosition(edgeCoord, 4);
428 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
429 interactionVolume.setEdgePosition(edgeCoord, 5);
430
431 break;
432 }
433 case 7:
434 {
435 DimVector edgeCoord(geometry.global(referenceElement.position(0, dim - 1)));
436 interactionVolume.setEdgePosition(edgeCoord, 1);
437 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
438 interactionVolume.setEdgePosition(edgeCoord, 4);
439 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
440 interactionVolume.setEdgePosition(edgeCoord, 3);
441
442 break;
443 }
444 }
445 }
446
447 // Choose the type of the hanging-node-interaction-volume depending on the number of stored fine
448 // elements (see dissertation M. Wolff, http://elib.uni-stuttgart.de/opus/volltexte/2013/8661/)
449 // and add missing elements and geometric information
450 switch (interactionVolume.getElementNumber())
451 {
452 // hanging-node interaction volume of type 5 or 7
453 case 2:
454 {
455 InteractionVolume hangingNodeVolume(problem_.gridView().grid());
456
457 std::vector<int> elemIdxOld;
458 for (int i = 0; i < InteractionVolume::subVolumeTotalNum; i++)
459 {
460 if (interactionVolume.hasSubVolumeElement(i))
461 elemIdxOld.push_back(i);
462 }
463
464 int zeroFaceIdx = IndexTranslator::getFaceIndexFromElements(elemIdxOld[0], elemIdxOld[1]);
465
466 std::vector<int> elemIdxNew(2);
467 elemIdxNew[0] = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elemIdxOld[0]);
468 elemIdxNew[1] = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elemIdxOld[1]);
469
470
471 hangingNodeVolume.setSubVolumeElement(interactionVolume.getSubVolumeElement(elemIdxOld[0]),
472 elemIdxNew[0]);
473 hangingNodeVolume.setSubVolumeElement(interactionVolume.getSubVolumeElement(elemIdxOld[1]),
474 elemIdxNew[1]);
475
476 for (int elem = 0; elem < InteractionVolume::subVolumeTotalNum; elem++)
477 {
478 for (int i = 0; i < 3; i++)
479 {
480 int faceIdxOld = IndexTranslator::getFaceIndexFromSubVolume(elem, i);
481 int faceIdxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, faceIdxOld);
482
483 for (int j = 0; j < 3; j++)
484 {
485 int elemNew = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elem);
486 int faceIdxNewTest = IndexTranslator::getFaceIndexFromSubVolume(elemNew, j);
487
488 if (faceIdxNew == faceIdxNewTest)
489 {
490 hangingNodeVolume.setNormal(interactionVolume.getNormal(elem, i),
491 elemNew, j);
492 hangingNodeVolume.setIndexOnElement(
493 interactionVolume.getIndexOnElement(elem, i), elemNew, j);
494 }
495 }
496 }
497 }
498
499 for (int i = 0; i < InteractionVolume::fluxFacesTotalNum; i++)
500 {
501 int idxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, i);
502 hangingNodeVolume.setFacePosition(interactionVolume.getFacePosition(i), idxNew);
503 }
504
505 for (int i = 0; i < InteractionVolume::fluxEdgesTotalNum; i++)
506 {
507 int idxNew = IndexTranslator::getNewEdgeIdxFromOldFaceIdxto0(zeroFaceIdx, i);
508 hangingNodeVolume.setEdgePosition(interactionVolume.getEdgePosition(i), idxNew);
509 }
510
511 interactionVolume = hangingNodeVolume;
512
513 interactionVolume.setHangingNodeType(InteractionVolume::twoSmallCells);
514
515 interactionVolume.setCenterPosition(centerPos);
516
517 auto element1 = interactionVolume.getSubVolumeElement(0);
518
519 for (const auto& intersection : intersections(problem_.gridView(), element1))
520 {
521 int idxInInside = intersection.indexInInside();
522
523 if (idxInInside == interactionVolume.getIndexOnElement(0, 2))
524 {
525 if (intersection.neighbor())
526 {
527 auto outside = intersection.outside();
528 if (element1.level() > outside.level())
529 {
530 interactionVolume.setSubVolumeElement(outside, 4);
531 interactionVolume.setSubVolumeElement(outside, 5);
532 }
533 }
534 }
535 else if (idxInInside == interactionVolume.getIndexOnElement(0, 1))
536 {
537 if (intersection.neighbor())
538 {
539 auto outside = intersection.outside();
540 if (element1.level() > outside.level())
541 {
542 interactionVolume.setSubVolumeElement(outside, 2);
543 interactionVolume.setSubVolumeElement(outside, 3);
544 }
545 }
546 }
547 }
548 auto element2 = interactionVolume.getSubVolumeElement(1);
549 auto element45 = interactionVolume.getSubVolumeElement(4);
550 auto element23 = interactionVolume.getSubVolumeElement(2);
551
552 for (const auto& intersection1 : intersections(problem_.gridView(), element45))
553 {
554 if (intersection1.neighbor())
555 {
556 auto element45Outside = intersection1.outside();
557
558 for (const auto& intersection2 : intersections(problem_.gridView(), element23))
559 {
560 if (intersection2.neighbor())
561 {
562 auto element23Outside = intersection2.outside();
563
564 if (element45Outside == element23Outside && element45Outside != element1
565 && element45Outside != element2)
566 {
567 interactionVolume.setSubVolumeElement(element45Outside, 6);
568 interactionVolume.setSubVolumeElement(element45Outside, 7);
569 DimVector normal = intersection2.centerUnitOuterNormal();
570 interactionVolume.setNormal(normal, 2, 2);
571 interactionVolume.setNormal(normal, 3, 2);
572 normal *= -1;
573 interactionVolume.setNormal(normal, 6, 0);
574 interactionVolume.setNormal(normal, 7, 0);
575
576 GlobalPosition globalPosFace = intersection2.geometry().center();
577 interactionVolume.setFacePosition(globalPosFace, 10);
578 interactionVolume.setFacePosition(globalPosFace, 11);
579 interactionVolume.setEdgePosition(centerPos, 4);
580
581 Scalar faceArea = intersection2.geometry().volume()/4.0;
582
583 interactionVolume.setFaceArea(faceArea, 10);
584 interactionVolume.setFaceArea(faceArea, 11);
585
586 normal = intersection1.centerUnitOuterNormal();
587 interactionVolume.setNormal(normal, 4, 2);
588 interactionVolume.setNormal(normal, 5, 1);
589 normal *= -1;
590 interactionVolume.setNormal(normal, 6, 1);
591 interactionVolume.setNormal(normal, 7, 2);
592
593 globalPosFace = intersection1.geometry().center();
594 interactionVolume.setFacePosition(globalPosFace, 5);
595 interactionVolume.setFacePosition(globalPosFace, 7);
596 interactionVolume.setEdgePosition(centerPos, 1);
597
598 faceArea = intersection1.geometry().volume()/4.0;
599
600 interactionVolume.setFaceArea(faceArea, 5);
601 interactionVolume.setFaceArea(faceArea, 7);
602 }
603 }
604 }
605 }
606 }
607
608 DimVector edgeCoord1(interactionVolume.getEdgePosition(0));
609 DimVector edgeCoord3(interactionVolume.getEdgePosition(2));
610 DimVector edgeCoord4(interactionVolume.getEdgePosition(3));
611 DimVector edgeCoord6(interactionVolume.getEdgePosition(5));
612
613 DimVector crossProductVector1(0);
614 DimVector crossProductVector2(0);
615
616 GlobalPosition globalPosFace = interactionVolume.getFacePosition(0);
617 crossProductVector1 = centerPos - globalPosFace;
618 crossProductVector2 = edgeCoord1 - edgeCoord3;
619 Scalar faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
620 interactionVolume.setFaceArea(faceArea, 0);
621
622 globalPosFace = interactionVolume.getFacePosition(1);
623 crossProductVector1 = centerPos - globalPosFace;
624 crossProductVector2 = edgeCoord1 - edgeCoord4;
625 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
626 interactionVolume.setFaceArea(faceArea, 1);
627
628 globalPosFace = interactionVolume.getFacePosition(3);
629 crossProductVector1 = centerPos - globalPosFace;
630 crossProductVector2 = edgeCoord1 - edgeCoord6;
631 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
632 interactionVolume.setFaceArea(faceArea, 3);
633
634 globalPosFace = interactionVolume.getFacePosition(8);
635 crossProductVector1 = centerPos - globalPosFace;
636 crossProductVector2 = edgeCoord3 - edgeCoord6;
637 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
638 interactionVolume.setFaceArea(faceArea, 8);
639
640 globalPosFace = interactionVolume.getFacePosition(9);
641 crossProductVector1 = centerPos - globalPosFace;
642 crossProductVector2 = edgeCoord3 - edgeCoord4;
643 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
644 interactionVolume.setFaceArea(faceArea, 9);
645
646 break;
647 }
648 //hanging-node interaction volume of type 1, 3 or 4
649 case 4:
650 {
651 InteractionVolume hangingNodeVolume(problem_.gridView().grid());
652
653 std::vector<int> elemIdxOld;
654 for (int i = 0; i < InteractionVolume::subVolumeTotalNum; i++)
655 {
656 if (interactionVolume.hasSubVolumeElement(i))
657 elemIdxOld.push_back(i);
658 }
659
660 std::set<int> zeroFaceIdxVec;
661 for (int i = 0; i < 4; i++)
662 {
663 for (int j = 0; j < 4; j++)
664 {
665 int fIdx = IndexTranslator::getFaceIndexFromElements(elemIdxOld[i], elemIdxOld[j]);
666 if (fIdx >= 0)
667 zeroFaceIdxVec.insert(fIdx);
668 }
669 }
670
671 std::vector<int> elemIdxNew(4);
672 int zeroFaceIdx = 0;
673
674 if (zeroFaceIdxVec.size() == 2)
675 {
676 hangingNodeVolume.setHangingNodeType(InteractionVolume::fourSmallCellsDiag);
677
678 if (zeroFaceIdxVec.find(0) == zeroFaceIdxVec.end())
679 zeroFaceIdx = *zeroFaceIdxVec.begin();
680
681 for (int i = 0; i < 4; i++)
682 {
683 elemIdxNew[i] = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elemIdxOld[i]);
684 hangingNodeVolume.setSubVolumeElement(interactionVolume.getSubVolumeElement(elemIdxOld[i]),
685 elemIdxNew[i]);
686 }
687 }
688 else if (zeroFaceIdxVec.size() == 4)
689 {
690 std::set<int>::iterator it = zeroFaceIdxVec.begin();
691 for (; it != zeroFaceIdxVec.end(); ++it)
692 {
693 zeroFaceIdx = *it;
694
695 bool isFace = true;
696 for (int i = 0; i < 4; i++)
697 {
698 elemIdxNew[i] = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx,
699 elemIdxOld[i]);
700 if (elemIdxNew[i] == 4 || elemIdxNew[i] == 5 || elemIdxNew[i] == 6 || elemIdxNew[i] == 7)
701 {
702 isFace = false;
703 break;
704 }
705 }
706 if (isFace)
707 {
708 for (int i = 0; i < 4; i++)
709 {
710 hangingNodeVolume.setSubVolumeElement(interactionVolume.getSubVolumeElement(elemIdxOld[i]), elemIdxNew[i]);
711 }
712
713 break;
714 }
715 }
716 }
717 else
718 {
719 std::set<int>::iterator it = zeroFaceIdxVec.begin();
720 std::cout<<"zeroFaceIdxVec = ";
721 for(;it != zeroFaceIdxVec.end(); ++it)
722 {
723 std::cout<<" "<<*it;
724 }
725 std::cout<<"\n";
726 }
727
728 for (int elem = 0; elem < InteractionVolume::subVolumeTotalNum; elem++)
729 {
730 for (int i = 0; i < 3; i++)
731 {
732 int faceIdxOld = IndexTranslator::getFaceIndexFromSubVolume(elem, i);
733 int faceIdxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, faceIdxOld);
734
735 for (int j = 0; j < 3; j++)
736 {
737 int elemNew = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elem);
738 int faceIdxNewTest = IndexTranslator::getFaceIndexFromSubVolume(elemNew, j);
739
740 if (faceIdxNew == faceIdxNewTest)
741 {
742 hangingNodeVolume.setNormal(interactionVolume.getNormal(elem, i),
743 elemNew, j);
744 hangingNodeVolume.setIndexOnElement(interactionVolume.getIndexOnElement(elem, i), elemNew, j);
745 }
746 }
747 }
748 }
749
750 for (int i = 0; i < InteractionVolume::fluxFacesTotalNum; i++)
751 {
752 int idxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, i);
753 hangingNodeVolume.setFacePosition(interactionVolume.getFacePosition(i), idxNew);
754 }
755
756 for (int i = 0; i < InteractionVolume::fluxEdgesTotalNum; i++)
757 {
758 int idxNew = IndexTranslator::getNewEdgeIdxFromOldFaceIdxto0(zeroFaceIdx, i);
759 hangingNodeVolume.setEdgePosition(interactionVolume.getEdgePosition(i), idxNew);
760 }
761
762 interactionVolume = hangingNodeVolume;
763
764 interactionVolume.setCenterPosition(centerPos);
765
766 if (zeroFaceIdxVec.size() == 4)
767 {
768 auto element1 = interactionVolume.getSubVolumeElement(0);
769 auto element4 = interactionVolume.getSubVolumeElement(3);
770
771 auto outside1 = element1;
772 auto outside4 = element4;
773
774 for (const auto& intersection1 : intersections(problem_.gridView(), element1))
775 {
776 if (intersection1.neighbor())
777 {
778 if (intersection1.indexInInside() == interactionVolume.getIndexOnElement(0, 2))
779 {
780 outside1 = intersection1.outside();
781 break;
782 }
783 }
784 }
785 for (const auto& intersection4 : intersections(problem_.gridView(), element4))
786 {
787 if (intersection4.neighbor())
788 {
789 if (intersection4.indexInInside() == interactionVolume.getIndexOnElement(3, 2))
790 {
791 outside4 = intersection4.outside();
792 break;
793 }
794 }
795 }
796 if (outside1 != outside4)
797 {
798 interactionVolume.setHangingNodeType(InteractionVolume::fourSmallCellsEdge);
799 }
800 else if (outside1 == outside4)
801 {
802 interactionVolume.setHangingNodeType(InteractionVolume::fourSmallCellsFace);
803 }
804 }
805
806 if (interactionVolume.getHangingNodeType() == InteractionVolume::fourSmallCellsFace)
807 {
808 auto element = interactionVolume.getSubVolumeElement(0);
809
810 DimVector edgeCoord1(interactionVolume.getEdgePosition(0));
811 DimVector edgeCoord3(interactionVolume.getEdgePosition(2));
812 DimVector edgeCoord4(interactionVolume.getEdgePosition(3));
813 DimVector edgeCoord5(interactionVolume.getEdgePosition(4));
814 DimVector edgeCoord6(interactionVolume.getEdgePosition(5));
815
816 DimVector crossProductVector1(0);
817 DimVector crossProductVector2(0);
818
819 GlobalPosition globalPosFace = interactionVolume.getFacePosition(0);
820 crossProductVector1 = centerPos - globalPosFace;
821 crossProductVector2 = edgeCoord1 - edgeCoord3;
822 Scalar faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
823 interactionVolume.setFaceArea(faceArea, 0);
824
825 globalPosFace = interactionVolume.getFacePosition(1);
826 crossProductVector1 = centerPos - globalPosFace;
827 crossProductVector2 = edgeCoord1 - edgeCoord4;
828 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
829 interactionVolume.setFaceArea(faceArea, 1);
830
831 globalPosFace = interactionVolume.getFacePosition(2);
832 crossProductVector1 = centerPos - globalPosFace;
833 crossProductVector2 = edgeCoord1 - edgeCoord5;
834 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
835 interactionVolume.setFaceArea(faceArea, 2);
836
837 globalPosFace = interactionVolume.getFacePosition(3);
838 crossProductVector1 = centerPos - globalPosFace;
839 crossProductVector2 = edgeCoord1 - edgeCoord6;
840 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
841 interactionVolume.setFaceArea(faceArea, 3);
842
843 globalPosFace = interactionVolume.getFacePosition(8);
844 crossProductVector1 = centerPos - globalPosFace;
845 crossProductVector2 = edgeCoord3 - edgeCoord6;
846 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
847 interactionVolume.setFaceArea(faceArea, 8);
848
849 globalPosFace = interactionVolume.getFacePosition(9);
850 crossProductVector1 = centerPos - globalPosFace;
851 crossProductVector2 = edgeCoord3 - edgeCoord4;
852 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
853 interactionVolume.setFaceArea(faceArea, 9);
854
855 globalPosFace = interactionVolume.getFacePosition(10);
856 crossProductVector1 = centerPos - globalPosFace;
857 crossProductVector2 = edgeCoord4 - edgeCoord5;
858 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
859 interactionVolume.setFaceArea(faceArea, 10);
860
861 globalPosFace = interactionVolume.getFacePosition(11);
862 crossProductVector1 = centerPos - globalPosFace;
863 crossProductVector2 = edgeCoord5 - edgeCoord6;
864 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
865 interactionVolume.setFaceArea(faceArea, 11);
866
867 for (const auto& intersection : intersections(problem_.gridView(), element))
868 {
869 if (intersection.indexInInside() == interactionVolume.getIndexOnElement(0, 2))
870 {
871 if (intersection.neighbor())
872 {
873 auto outside = intersection.outside();
874 if (element.level() > outside.level())
875 {
876 interactionVolume.setSubVolumeElement(outside, 4);
877 interactionVolume.setSubVolumeElement(outside, 5);
878 interactionVolume.setSubVolumeElement(outside, 6);
879 interactionVolume.setSubVolumeElement(outside, 7);
880
881 break;
882 }
883 }
884 }
885
886 }
887 }
888 else if (interactionVolume.getHangingNodeType() == InteractionVolume::fourSmallCellsEdge)
889 {
890 auto element1 = interactionVolume.getSubVolumeElement(0);
891 auto element2 = interactionVolume.getSubVolumeElement(1);
892 auto element3 = interactionVolume.getSubVolumeElement(2);
893 auto element4 = interactionVolume.getSubVolumeElement(3);
894
895 for (const auto& intersection1 : intersections(problem_.gridView(), element1))
896 {
897 if (intersection1.neighbor())
898 {
899 if (intersection1.indexInInside() == interactionVolume.getIndexOnElement(0, 2))
900 {
901 auto outside = intersection1.outside();
902 if (element1.level() > outside.level())
903 {
904 interactionVolume.setSubVolumeElement(outside, 4);
905 }
906 }
907 }
908 }
909 for (const auto& intersection2 : intersections(problem_.gridView(), element2))
910 {
911 if (intersection2.neighbor())
912 {
913 if (intersection2.indexInInside() == interactionVolume.getIndexOnElement(1, 2))
914 {
915 auto outside = intersection2.outside();
916 if (element2.level() > outside.level())
917 {
918 interactionVolume.setSubVolumeElement(outside, 5);
919
920 break;
921 }
922 }
923 }
924 }
925 for (const auto& intersection3 : intersections(problem_.gridView(), element3))
926 {
927 if (intersection3.neighbor())
928 {
929 if (intersection3.indexInInside() == interactionVolume.getIndexOnElement(2, 2))
930 {
931 auto outside = intersection3.outside();
932 if (element3.level() > outside.level())
933 {
934 interactionVolume.setSubVolumeElement(outside, 6);
935 break;
936 }
937 }
938 }
939 }
940 for (const auto& intersection4 : intersections(problem_.gridView(), element4))
941 {
942 if (intersection4.neighbor())
943 {
944 if (intersection4.indexInInside() == interactionVolume.getIndexOnElement(3, 2))
945 {
946 auto outside = intersection4.outside();
947 if (element4.level() > outside.level())
948 {
949 interactionVolume.setSubVolumeElement(outside, 7);
950
951 break;
952 }
953 }
954 }
955 }
956
957 DimVector edgeCoord1(interactionVolume.getEdgePosition(0));
958 DimVector edgeCoord3(interactionVolume.getEdgePosition(2));
959 DimVector edgeCoord4(interactionVolume.getEdgePosition(3));
960 DimVector edgeCoord5(interactionVolume.getEdgePosition(4));
961 DimVector edgeCoord6(interactionVolume.getEdgePosition(5));
962
963 DimVector crossProductVector1(0);
964 DimVector crossProductVector2(0);
965
966 GlobalPosition globalPosFace = interactionVolume.getFacePosition(0);
967 crossProductVector1 = centerPos - globalPosFace;
968 crossProductVector2 = edgeCoord1 - edgeCoord3;
969 Scalar faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
970 interactionVolume.setFaceArea(faceArea, 0);
971
972 globalPosFace = interactionVolume.getFacePosition(1);
973 crossProductVector1 = centerPos - globalPosFace;
974 crossProductVector2 = edgeCoord1 - edgeCoord4;
975 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
976 interactionVolume.setFaceArea(faceArea, 1);
977
978 globalPosFace = interactionVolume.getFacePosition(2);
979 crossProductVector1 = centerPos - globalPosFace;
980 crossProductVector2 = edgeCoord1 - edgeCoord5;
981 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
982 interactionVolume.setFaceArea(faceArea, 2);
983
984 globalPosFace = interactionVolume.getFacePosition(3);
985 crossProductVector1 = centerPos - globalPosFace;
986 crossProductVector2 = edgeCoord1 - edgeCoord6;
987 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
988 interactionVolume.setFaceArea(faceArea, 3);
989
990 globalPosFace = interactionVolume.getFacePosition(8);
991 crossProductVector1 = centerPos - globalPosFace;
992 crossProductVector2 = edgeCoord3 - edgeCoord6;
993 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
994 interactionVolume.setFaceArea(faceArea, 8);
995
996 globalPosFace = interactionVolume.getFacePosition(9);
997 crossProductVector1 = centerPos - globalPosFace;
998 crossProductVector2 = edgeCoord3 - edgeCoord4;
999 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1000 interactionVolume.setFaceArea(faceArea, 9);
1001
1002 globalPosFace = interactionVolume.getFacePosition(10);
1003 crossProductVector1 = centerPos - globalPosFace;
1004 crossProductVector2 = edgeCoord4 - edgeCoord5;
1005 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1006 interactionVolume.setFaceArea(faceArea, 10);
1007
1008 globalPosFace = interactionVolume.getFacePosition(11);
1009 crossProductVector1 = centerPos - globalPosFace;
1010 crossProductVector2 = edgeCoord5 - edgeCoord6;
1011 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1012 interactionVolume.setFaceArea(faceArea, 11);
1013
1014 auto element5 = interactionVolume.getSubVolumeElement(4);
1015 auto element6 = interactionVolume.getSubVolumeElement(5);
1016 auto element7 = interactionVolume.getSubVolumeElement(6);
1017 auto element8 = interactionVolume.getSubVolumeElement(7);
1018
1019 if (element5 == element6)
1020 {
1021 interactionVolume.setFacePosition(element5.geometry().center(), 4);
1022 interactionVolume.setFacePosition(element7.geometry().center(), 6);
1023
1024 for (const auto& intersection : intersections(problem_.gridView(), element5))
1025 {
1026 if (intersection.neighbor())
1027 {
1028 auto outside = intersection.outside();
1029
1030 if (outside == element7 || outside == element8)
1031 {
1032 int indexInInside = intersection.indexInInside();
1033 interactionVolume.setIndexOnElement(indexInInside, 4, 2);
1034 interactionVolume.setIndexOnElement(indexInInside, 5, 1);
1035 DimVector normal = intersection.centerUnitOuterNormal();
1036 interactionVolume.setNormal(normal, 4, 2);
1037 interactionVolume.setNormal(normal, 5, 1);
1038 globalPosFace = intersection.geometry().center();
1039 interactionVolume.setFacePosition(globalPosFace, 5);
1040 interactionVolume.setFacePosition(globalPosFace, 7);
1041 int indexInOutside = intersection.indexInOutside();
1042 interactionVolume.setIndexOnElement(indexInOutside, 6, 1);
1043 interactionVolume.setIndexOnElement(indexInOutside, 7, 2);
1044 normal *= -1;
1045 interactionVolume.setNormal(normal, 6, 1);
1046 interactionVolume.setNormal(normal, 7, 2);
1047
1048 break;
1049 }
1050 }
1051 }
1052
1053 DimVector edgeCoord2(interactionVolume.getFacePosition(7));
1054
1055 globalPosFace = interactionVolume.getFacePosition(5);
1056 crossProductVector1 = centerPos - globalPosFace;
1057 crossProductVector2 = edgeCoord2 - edgeCoord4;
1058 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1059 interactionVolume.setFaceArea(faceArea, 5);
1060 // interactionVolume.setFaceArea(0.0, 5);
1061 globalPosFace = interactionVolume.getFacePosition(7);
1062 crossProductVector1 = centerPos - globalPosFace;
1063 crossProductVector2 = edgeCoord2 - edgeCoord6;
1064 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1065 interactionVolume.setFaceArea(faceArea, 7);
1066 // interactionVolume.setFaceArea(0.0, 7);
1067 }
1068 else if (element5 == element7)
1069 {
1070 interactionVolume.setFacePosition(element6.geometry().center(), 5);
1071 interactionVolume.setFacePosition(element5.geometry().center(), 7);
1072 interactionVolume.setFacePosition(globalPosFace, 6);
1073
1074 for (const auto& intersection : intersections(problem_.gridView(), element5))
1075 {
1076 if (intersection.neighbor())
1077 {
1078 auto outside = intersection.outside();
1079
1080 if (outside == element6 || outside == element8)
1081 {
1082 int indexInInside = intersection.indexInInside();
1083 interactionVolume.setIndexOnElement(indexInInside, 4, 1);
1084 interactionVolume.setIndexOnElement(indexInInside, 6, 2);
1085 DimVector normal = intersection.centerUnitOuterNormal();
1086 interactionVolume.setNormal(normal, 4, 1);
1087 interactionVolume.setNormal(normal, 6, 2);
1088 globalPosFace = intersection.geometry().center();
1089 interactionVolume.setFacePosition(globalPosFace, 4);
1090 interactionVolume.setFacePosition(globalPosFace, 6);
1091 int indexInOutside = intersection.indexInOutside();
1092 interactionVolume.setIndexOnElement(indexInOutside, 5, 2);
1093 interactionVolume.setIndexOnElement(indexInOutside, 7, 1);
1094 normal *= -1;
1095 interactionVolume.setNormal(normal, 5, 2);
1096 interactionVolume.setNormal(normal, 7, 1);
1097
1098 break;
1099 }
1100 }
1101 }
1102 DimVector edgeCoord2(interactionVolume.getFacePosition(4));
1103
1104 globalPosFace = interactionVolume.getFacePosition(4);
1105 crossProductVector1 = centerPos - globalPosFace;
1106 crossProductVector2 = edgeCoord2 - edgeCoord3;
1107 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1108 interactionVolume.setFaceArea(faceArea, 4);
1109
1110 globalPosFace = interactionVolume.getFacePosition(6);
1111 crossProductVector1 = centerPos - globalPosFace;
1112 crossProductVector2 = edgeCoord2 - edgeCoord5;
1113 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1114 interactionVolume.setFaceArea(faceArea, 6);
1115 // interactionVolume.setFaceArea(0.0, 4);
1116 // interactionVolume.setFaceArea(0.0, 6);
1117 }
1118
1119 const ElementGeometry& geometry = element5.geometry();
1120
1121 const auto referenceElement = ReferenceElements::general(geometry.type());
1122
1123 int oldSubVolumElemIdx = IndexTranslator::getOldElemIdxFromNewFaceIdxto0(zeroFaceIdx, 4);
1124 int oldEdgeIdx = IndexTranslator::getOldEdgeIdxFromNewFaceIdxto0(zeroFaceIdx, 1);
1125
1126 DimVector edgeCoord(0.0);
1127 switch (oldSubVolumElemIdx)
1128 {
1129 case 0:
1130 {
1131 switch (oldEdgeIdx)
1132 {
1133 case 2:
1134 edgeCoord = geometry.global(referenceElement.position(9, dim - 1));
1135 break;
1136 case 0:
1137 edgeCoord = geometry.global(referenceElement.position(3, dim - 1));
1138 break;
1139 case 5:
1140 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
1141 break;
1142 }
1143
1144 break;
1145 }
1146 case 1:
1147 {
1148 switch (oldEdgeIdx)
1149 {
1150 case 0:
1151 edgeCoord = geometry.global(referenceElement.position(2, dim - 1));
1152 break;
1153 case 2:
1154 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
1155 break;
1156 case 3:
1157 edgeCoord = geometry.global(referenceElement.position(11, dim - 1));
1158 break;
1159 }
1160
1161 break;
1162 }
1163 case 2:
1164 {
1165 switch (oldEdgeIdx)
1166 {
1167 case 0:
1168 edgeCoord = geometry.global(referenceElement.position(1, dim - 1));
1169 break;
1170 case 4:
1171 edgeCoord = geometry.global(referenceElement.position(9, dim - 1));
1172 break;
1173 case 5:
1174 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
1175 break;
1176 }
1177
1178 break;
1179 }
1180 case 3:
1181 {
1182 switch (oldEdgeIdx)
1183 {
1184 case 0:
1185 edgeCoord = geometry.global(referenceElement.position(0, dim - 1));
1186 break;
1187 case 4:
1188 edgeCoord = geometry.global(referenceElement.position(8, dim - 1));
1189 break;
1190 case 3:
1191 edgeCoord = geometry.global(referenceElement.position(10, dim - 1));
1192 break;
1193 }
1194
1195 break;
1196 }
1197 case 4:
1198 {
1199 switch (oldEdgeIdx)
1200 {
1201 case 1:
1202 edgeCoord = geometry.global(referenceElement.position(3, dim - 1));
1203 break;
1204 case 2:
1205 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
1206 break;
1207 case 5:
1208 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
1209 break;
1210 }
1211
1212 break;
1213 }
1214 case 5:
1215 {
1216 switch (oldEdgeIdx)
1217 {
1218 case 1:
1219 edgeCoord = geometry.global(referenceElement.position(2, dim - 1));
1220 break;
1221 case 2:
1222 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
1223 break;
1224 case 3:
1225 edgeCoord = geometry.global(referenceElement.position(7, dim - 1));
1226 break;
1227 }
1228
1229 break;
1230 }
1231 case 6:
1232 {
1233 switch (oldEdgeIdx)
1234 {
1235 case 1:
1236 edgeCoord = geometry.global(referenceElement.position(1, dim - 1));
1237 break;
1238 case 4:
1239 edgeCoord = geometry.global(referenceElement.position(5, dim - 1));
1240 break;
1241 case 5:
1242 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
1243 break;
1244 }
1245
1246 break;
1247 }
1248 case 7:
1249 {
1250 switch (oldEdgeIdx)
1251 {
1252 case 1:
1253 edgeCoord = geometry.global(referenceElement.position(0, dim - 1));
1254 break;
1255 case 4:
1256 edgeCoord = geometry.global(referenceElement.position(4, dim - 1));
1257 break;
1258 case 3:
1259 edgeCoord = geometry.global(referenceElement.position(6, dim - 1));
1260 break;
1261 }
1262
1263 break;
1264 }
1265 }
1266
1267 interactionVolume.setEdgePosition(edgeCoord, 1);
1268 }
1269 else if (interactionVolume.getHangingNodeType() == InteractionVolume::fourSmallCellsDiag)
1270 {
1271 DimVector edgeCoord1(interactionVolume.getEdgePosition(0));
1272 DimVector edgeCoord2(interactionVolume.getEdgePosition(1));
1273 DimVector edgeCoord3(interactionVolume.getEdgePosition(2));
1274 DimVector edgeCoord4(interactionVolume.getEdgePosition(3));
1275 DimVector edgeCoord5(interactionVolume.getEdgePosition(4));
1276 DimVector edgeCoord6(interactionVolume.getEdgePosition(5));
1277
1278 DimVector crossProductVector1(0);
1279 DimVector crossProductVector2(0);
1280
1281 GlobalPosition globalPosFace = interactionVolume.getFacePosition(0);
1282 crossProductVector1 = centerPos - globalPosFace;
1283 crossProductVector2 = edgeCoord1 - edgeCoord3;
1284 Scalar faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1285 interactionVolume.setFaceArea(faceArea, 0);
1286
1287 globalPosFace = interactionVolume.getFacePosition(1);
1288 crossProductVector1 = centerPos - globalPosFace;
1289 crossProductVector2 = edgeCoord1 - edgeCoord4;
1290 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1291 interactionVolume.setFaceArea(faceArea, 1);
1292
1293 globalPosFace = interactionVolume.getFacePosition(3);
1294 crossProductVector1 = centerPos - globalPosFace;
1295 crossProductVector2 = edgeCoord1 - edgeCoord6;
1296 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1297 interactionVolume.setFaceArea(faceArea, 3);
1298
1299 globalPosFace = interactionVolume.getFacePosition(5);
1300 crossProductVector1 = centerPos - globalPosFace;
1301 crossProductVector2 = edgeCoord2 - edgeCoord4;
1302 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1303 interactionVolume.setFaceArea(faceArea, 5);
1304
1305 globalPosFace = interactionVolume.getFacePosition(6);
1306 crossProductVector1 = centerPos - globalPosFace;
1307 crossProductVector2 = edgeCoord2 - edgeCoord5;
1308 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1309 interactionVolume.setFaceArea(faceArea, 6);
1310
1311 globalPosFace = interactionVolume.getFacePosition(7);
1312 crossProductVector1 = centerPos - globalPosFace;
1313 crossProductVector2 = edgeCoord2 - edgeCoord6;
1314 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1315 interactionVolume.setFaceArea(faceArea, 7);
1316
1317 globalPosFace = interactionVolume.getFacePosition(8);
1318 crossProductVector1 = centerPos - globalPosFace;
1319 crossProductVector2 = edgeCoord3 - edgeCoord6;
1320 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1321 interactionVolume.setFaceArea(faceArea, 8);
1322
1323 globalPosFace = interactionVolume.getFacePosition(9);
1324 crossProductVector1 = centerPos - globalPosFace;
1325 crossProductVector2 = edgeCoord3 - edgeCoord4;
1326 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1327 interactionVolume.setFaceArea(faceArea, 9);
1328
1329 globalPosFace = interactionVolume.getFacePosition(10);
1330 crossProductVector1 = centerPos - globalPosFace;
1331 crossProductVector2 = edgeCoord4 - edgeCoord5;
1332 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1333 interactionVolume.setFaceArea(faceArea, 10);
1334
1335 globalPosFace = interactionVolume.getFacePosition(11);
1336 crossProductVector1 = centerPos - globalPosFace;
1337 crossProductVector2 = edgeCoord5 - edgeCoord6;
1338 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1339 interactionVolume.setFaceArea(faceArea, 11);
1340
1341 auto element1 = interactionVolume.getSubVolumeElement(0);
1342
1343 bool hasFaceOne = false;
1344 bool hasFaceTwo = false;
1345 for (const auto& intersection : intersections(problem_.gridView(), element1))
1346 {
1347 if (intersection.indexInInside() == interactionVolume.getIndexOnElement(0, 1))
1348 {
1349 if (intersection.neighbor())
1350 {
1351 auto outside = intersection.outside();
1352 if (element1.level() > outside.level())
1353 {
1354 interactionVolume.setSubVolumeElement(outside, 2);
1355 interactionVolume.setSubVolumeElement(outside, 3);
1356
1357 hasFaceOne = true;
1358 if (hasFaceTwo)
1359 break;
1360 }
1361 }
1362 }
1363 if (intersection.indexInInside() == interactionVolume.getIndexOnElement(0, 2))
1364 {
1365 if (intersection.neighbor())
1366 {
1367 auto outside = intersection.outside();
1368 if (element1.level() > outside.level())
1369 {
1370 interactionVolume.setSubVolumeElement(outside, 4);
1371 interactionVolume.setSubVolumeElement(outside, 5);
1372
1373 hasFaceTwo = true;
1374 if (hasFaceOne)
1375 break;
1376 }
1377 }
1378 }
1379 }
1380 }
1381
1382 break;
1383 }
1384 //hanging-node interaction volume of type 2 or 6
1385 case 6:
1386 {
1387 InteractionVolume hangingNodeVolume(problem_.gridView().grid());
1388
1389 std::vector<int> elemIdxOld;
1390 for (int i = 0; i < InteractionVolume::subVolumeTotalNum; i++)
1391 {
1392 if (interactionVolume.hasSubVolumeElement(i))
1393 elemIdxOld.push_back(i);
1394 }
1395
1396 std::set<int> zeroFaceIdxVec;
1397 for (int i = 0; i < 6; i++)
1398 {
1399 for (int j = 0; j < 6; j++)
1400 {
1401 int fIdx = IndexTranslator::getFaceIndexFromElements(elemIdxOld[i], elemIdxOld[j]);
1402 if (fIdx >= 0)
1403 zeroFaceIdxVec.insert(fIdx);
1404 }
1405 }
1406
1407 std::vector<int> elemIdxNew(6);
1408 int zeroFaceIdx = 0;
1409
1410 bool isFace = true;
1411 std::set<int>::iterator it = zeroFaceIdxVec.begin();
1412 for (; it != zeroFaceIdxVec.end(); ++it)
1413 {
1414 zeroFaceIdx = *it;
1415 isFace = true;
1416 // bool isFace = true;
1417 for (int i = 0; i < 6; i++)
1418 {
1419 elemIdxNew[i] = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elemIdxOld[i]);
1420 if (elemIdxNew[i] == 6 || elemIdxNew[i] == 7)
1421 {
1422 isFace = false;
1423 break;
1424 }
1425 }
1426 if (isFace)
1427 {
1428 for (int i = 0; i < 6; i++)
1429 {
1430 hangingNodeVolume.setSubVolumeElement(interactionVolume.getSubVolumeElement(elemIdxOld[i]),
1431 elemIdxNew[i]);
1432 }
1433
1434 break;
1435 }
1436 }
1437
1438 for (int elem = 0; elem < InteractionVolume::subVolumeTotalNum; elem++)
1439 {
1440 for (int i = 0; i < 3; i++)
1441 {
1442 int faceIdxOld = IndexTranslator::getFaceIndexFromSubVolume(elem, i);
1443 int faceIdxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, faceIdxOld);
1444
1445 for (int j = 0; j < 3; j++)
1446 {
1447 int elemNew = IndexTranslator::getNewElemIdxFromOldFaceIdxto0(zeroFaceIdx, elem);
1448 int faceIdxNewTest = IndexTranslator::getFaceIndexFromSubVolume(elemNew, j);
1449
1450 if (faceIdxNew == faceIdxNewTest)
1451 {
1452 hangingNodeVolume.setNormal(interactionVolume.getNormal(elem, i),
1453 elemNew, j);
1454 hangingNodeVolume.setIndexOnElement(
1455 interactionVolume.getIndexOnElement(elem, i), elemNew, j);
1456 }
1457 }
1458 }
1459 }
1460
1461 for (int i = 0; i < InteractionVolume::fluxFacesTotalNum; i++)
1462 {
1463 int idxNew = IndexTranslator::getNewFaceIdxFromOldIdxto0(zeroFaceIdx, i);
1464 hangingNodeVolume.setFaceArea(interactionVolume.getFaceArea(i), idxNew);
1465 hangingNodeVolume.setFacePosition(interactionVolume.getFacePosition(i), idxNew);
1466 }
1467
1468 for (int i = 0; i < InteractionVolume::fluxEdgesTotalNum; i++)
1469 {
1470 int idxNew = IndexTranslator::getNewEdgeIdxFromOldFaceIdxto0(zeroFaceIdx, i);
1471 hangingNodeVolume.setEdgePosition(interactionVolume.getEdgePosition(i), idxNew);
1472 }
1473
1474 interactionVolume = hangingNodeVolume;
1475
1476 interactionVolume.setCenterPosition(centerPos);
1477
1478 interactionVolume.setHangingNodeType(InteractionVolume::sixSmallCells);
1479
1480 auto element3 = interactionVolume.getSubVolumeElement(2);
1481
1482 for (const auto& intersection : intersections(problem_.gridView(), element3))
1483 {
1484 if (intersection.indexInInside() == interactionVolume.getIndexOnElement(2, 2))
1485 {
1486 if (intersection.neighbor())
1487 {
1488 auto outside = intersection.outside();
1489 if (element3.level() > outside.level())
1490 {
1491 interactionVolume.setSubVolumeElement(outside, 6);
1492 interactionVolume.setSubVolumeElement(outside, 7);
1493
1494 break;
1495 }
1496 }
1497 }
1498 }
1499
1500 DimVector edgeCoord1(interactionVolume.getEdgePosition(0));
1501 DimVector edgeCoord2(interactionVolume.getEdgePosition(1));
1502 DimVector edgeCoord3(interactionVolume.getEdgePosition(2));
1503 DimVector edgeCoord4(interactionVolume.getEdgePosition(3));
1504 DimVector edgeCoord5(interactionVolume.getEdgePosition(4));
1505 DimVector edgeCoord6(interactionVolume.getEdgePosition(5));
1506
1507 DimVector crossProductVector1(0);
1508 DimVector crossProductVector2(0);
1509
1510 GlobalPosition globalPosFace = interactionVolume.getFacePosition(0);
1511 crossProductVector1 = centerPos - globalPosFace;
1512 crossProductVector2 = edgeCoord1 - edgeCoord3;
1513 Scalar faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1514 interactionVolume.setFaceArea(faceArea, 0);
1515
1516 globalPosFace = interactionVolume.getFacePosition(1);
1517 crossProductVector1 = centerPos - globalPosFace;
1518 crossProductVector2 = edgeCoord1 - edgeCoord4;
1519 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1520 interactionVolume.setFaceArea(faceArea, 1);
1521
1522 globalPosFace = interactionVolume.getFacePosition(2);
1523 crossProductVector1 = centerPos - globalPosFace;
1524 crossProductVector2 = edgeCoord1 - edgeCoord5;
1525 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1526 interactionVolume.setFaceArea(faceArea, 2);
1527
1528 globalPosFace = interactionVolume.getFacePosition(3);
1529 crossProductVector1 = centerPos - globalPosFace;
1530 crossProductVector2 = edgeCoord1 - edgeCoord6;
1531 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1532 interactionVolume.setFaceArea(faceArea, 3);
1533
1534 globalPosFace = interactionVolume.getFacePosition(4);
1535 crossProductVector1 = centerPos - globalPosFace;
1536 crossProductVector2 = edgeCoord2 - edgeCoord3;
1537 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1538 interactionVolume.setFaceArea(faceArea, 4);
1539
1540 globalPosFace = interactionVolume.getFacePosition(5);
1541 crossProductVector1 = centerPos - globalPosFace;
1542 crossProductVector2 = edgeCoord2 - edgeCoord4;
1543 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1544 interactionVolume.setFaceArea(faceArea, 5);
1545
1546 globalPosFace = interactionVolume.getFacePosition(7);
1547 crossProductVector1 = centerPos - globalPosFace;
1548 crossProductVector2 = edgeCoord2 - edgeCoord6;
1549 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1550 interactionVolume.setFaceArea(faceArea, 7);
1551
1552 globalPosFace = interactionVolume.getFacePosition(8);
1553 crossProductVector1 = centerPos - globalPosFace;
1554 crossProductVector2 = edgeCoord3 - edgeCoord6;
1555 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1556 interactionVolume.setFaceArea(faceArea, 8);
1557
1558 globalPosFace = interactionVolume.getFacePosition(9);
1559 crossProductVector1 = centerPos - globalPosFace;
1560 crossProductVector2 = edgeCoord3 - edgeCoord4;
1561 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1562 interactionVolume.setFaceArea(faceArea, 9);
1563
1564 globalPosFace = interactionVolume.getFacePosition(10);
1565 crossProductVector1 = centerPos - globalPosFace;
1566 crossProductVector2 = edgeCoord4 - edgeCoord5;
1567 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1568 interactionVolume.setFaceArea(faceArea, 10);
1569
1570 globalPosFace = interactionVolume.getFacePosition(11);
1571 crossProductVector1 = centerPos - globalPosFace;
1572 crossProductVector2 = edgeCoord5 - edgeCoord6;
1573 faceArea = crossProduct(crossProductVector1, crossProductVector2).two_norm()/2.0;
1574 interactionVolume.setFaceArea(faceArea, 11);
1575 break;
1576 }
1577 default:
1578 {
1579 interactionVolume.printInteractionVolumeInfo();
1580 DUNE_THROW(Dune::NotImplemented, "Hanging node shape not implemented");
1581 break;
1582 }
1583 }
1584}
1585
1587template<class TypeTag>
1589{
1590 std::vector < std::vector<int> > elemVertMap(problem_.gridView().size(dim), std::vector<int>(8, -1));
1591
1592 //Add elements to the interaction volumes and store element-vertex map
1593 for (const auto& element : elements(problem_.gridView()))
1594 asImp_().storeSubVolumeElements(element, elemVertMap);
1595
1596 for (unsigned int i = 0; i < asImp_().interactionVolumes_.size(); i++)
1597 if (asImp_().interactionVolumes_[i].getElementNumber() == 0)
1598 asImp_().interactionVolumes_[i].printInteractionVolumeInfo();
1599
1600 // Store information related to DUNE intersections for all interaction volumes
1601 for (const auto& element : elements(problem_.gridView()))
1602 asImp_().storeIntersectionInfo(element, elemVertMap);
1603
1604 faceVertices_.clear();
1605 faceVertices_.resize(problem_.gridView().size(0));
1606
1607 // Complete storage of the interaction volumes using the previously stored information
1608 // about the orientation and relationship of the DUNE elements in the interaction volumes (see doc/docextra/3dmpfa)
1609 for (const auto& vertex : vertices(problem_.gridView()))
1610 {
1611 int vIdxGlobal = problem_.variables().index(vertex);
1612
1613 InteractionVolume& interactionVolume = asImp_().interactionVolumes_[vIdxGlobal];
1614
1615 if (interactionVolume.getElementNumber() == 8)
1616 {
1617 asImp_().storeInnerInteractionVolume(interactionVolume, vertex);
1618 }
1619 else if (interactionVolume.isBoundaryInteractionVolume())
1620 {
1621 asImp_().storeBoundaryInteractionVolume(interactionVolume, vertex);
1622 }
1623 //hanging node!
1624 else
1625 {
1626 storeHangingNodeInteractionVolume(interactionVolume, vertex);
1627 }
1628
1629 if (!interactionVolume.isBoundaryInteractionVolume())
1630 {
1631 auto element1 = interactionVolume.getSubVolumeElement(0);
1632 auto element2 = interactionVolume.getSubVolumeElement(1);
1633 auto element3 = interactionVolume.getSubVolumeElement(2);
1634 auto element4 = interactionVolume.getSubVolumeElement(3);
1635 auto element5 = interactionVolume.getSubVolumeElement(4);
1636 auto element6 = interactionVolume.getSubVolumeElement(5);
1637 auto element7 = interactionVolume.getSubVolumeElement(6);
1638 auto element8 = interactionVolume.getSubVolumeElement(7);
1639
1640 int globalIdx1 = problem_.variables().index(element1);
1641 int globalIdx2 = problem_.variables().index(element2);
1642 int globalIdx3 = problem_.variables().index(element3);
1643 int globalIdx4 = problem_.variables().index(element4);
1644 int globalIdx5 = problem_.variables().index(element5);
1645 int globalIdx6 = problem_.variables().index(element6);
1646 int globalIdx7 = problem_.variables().index(element7);
1647 int globalIdx8 = problem_.variables().index(element8);
1648
1649 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(0, 0), globalIdx1, interactionVolume.getIndexOnElement(0, 0));
1650 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(0, 1), globalIdx1, interactionVolume.getIndexOnElement(0, 1));
1651 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(0, 2), globalIdx1, interactionVolume.getIndexOnElement(0, 2));
1652 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(1, 0), globalIdx2, interactionVolume.getIndexOnElement(1, 0));
1653 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(1, 1), globalIdx2, interactionVolume.getIndexOnElement(1, 1));
1654 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(1, 2), globalIdx2, interactionVolume.getIndexOnElement(1, 2));
1655 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(2, 0), globalIdx3, interactionVolume.getIndexOnElement(2, 0));
1656 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(3, 1), globalIdx4, interactionVolume.getIndexOnElement(3, 1));
1657 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(4, 0), globalIdx5, interactionVolume.getIndexOnElement(4, 0));
1658 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(5, 0), globalIdx6, interactionVolume.getIndexOnElement(5, 0));
1659
1660 faceVertices_[globalIdx1][interactionVolume.getIndexOnElement(0, 0)].insert(vIdxGlobal);
1661 faceVertices_[globalIdx1][interactionVolume.getIndexOnElement(0, 1)].insert(vIdxGlobal);
1662 faceVertices_[globalIdx1][interactionVolume.getIndexOnElement(0, 2)].insert(vIdxGlobal);
1663 faceVertices_[globalIdx2][interactionVolume.getIndexOnElement(1, 0)].insert(vIdxGlobal);
1664 faceVertices_[globalIdx2][interactionVolume.getIndexOnElement(1, 1)].insert(vIdxGlobal);
1665 faceVertices_[globalIdx2][interactionVolume.getIndexOnElement(1, 2)].insert(vIdxGlobal);
1666 faceVertices_[globalIdx3][interactionVolume.getIndexOnElement(2, 0)].insert(vIdxGlobal);
1667 faceVertices_[globalIdx4][interactionVolume.getIndexOnElement(3, 1)].insert(vIdxGlobal);
1668 faceVertices_[globalIdx5][interactionVolume.getIndexOnElement(4, 0)].insert(vIdxGlobal);
1669 faceVertices_[globalIdx6][interactionVolume.getIndexOnElement(5, 0)].insert(vIdxGlobal);
1670
1671 if (interactionVolume.getHangingNodeType() != InteractionVolume::twoSmallCells)
1672 {
1673 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(2, 1), globalIdx3, interactionVolume.getIndexOnElement(2, 1));
1674 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(2, 2), globalIdx3, interactionVolume.getIndexOnElement(2, 2));
1675 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(3, 0), globalIdx4, interactionVolume.getIndexOnElement(3, 0));
1676 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(3, 2), globalIdx4, interactionVolume.getIndexOnElement(3, 2));
1677 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(4, 1), globalIdx5, interactionVolume.getIndexOnElement(4, 1));
1678 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(4, 2), globalIdx5, interactionVolume.getIndexOnElement(4, 2));
1679 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(5, 1), globalIdx6, interactionVolume.getIndexOnElement(5, 1));
1680 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(5, 2), globalIdx6, interactionVolume.getIndexOnElement(5, 2));
1681 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(6, 0), globalIdx7, interactionVolume.getIndexOnElement(6, 0));
1682 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(6, 1), globalIdx7, interactionVolume.getIndexOnElement(6, 1));
1683 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(6, 2), globalIdx7, interactionVolume.getIndexOnElement(6, 2));
1684 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(7, 0), globalIdx8, interactionVolume.getIndexOnElement(7, 0));
1685 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(7, 1), globalIdx8, interactionVolume.getIndexOnElement(7, 1));
1686 asImp_().addRealFluxFaceArea_(interactionVolume.getFaceArea(7, 2), globalIdx8, interactionVolume.getIndexOnElement(7, 2));
1687
1688 faceVertices_[globalIdx3][interactionVolume.getIndexOnElement(2, 1)].insert(vIdxGlobal);
1689 faceVertices_[globalIdx3][interactionVolume.getIndexOnElement(2, 2)].insert(vIdxGlobal);
1690 faceVertices_[globalIdx4][interactionVolume.getIndexOnElement(3, 0)].insert(vIdxGlobal);
1691 faceVertices_[globalIdx4][interactionVolume.getIndexOnElement(3, 2)].insert(vIdxGlobal);
1692 faceVertices_[globalIdx5][interactionVolume.getIndexOnElement(4, 1)].insert(vIdxGlobal);
1693 faceVertices_[globalIdx5][interactionVolume.getIndexOnElement(4, 2)].insert(vIdxGlobal);
1694 faceVertices_[globalIdx6][interactionVolume.getIndexOnElement(5, 1)].insert(vIdxGlobal);
1695 faceVertices_[globalIdx6][interactionVolume.getIndexOnElement(5, 2)].insert(vIdxGlobal);
1696 faceVertices_[globalIdx7][interactionVolume.getIndexOnElement(6, 0)].insert(vIdxGlobal);
1697 faceVertices_[globalIdx7][interactionVolume.getIndexOnElement(6, 1)].insert(vIdxGlobal);
1698 faceVertices_[globalIdx7][interactionVolume.getIndexOnElement(6, 2)].insert(vIdxGlobal);
1699 faceVertices_[globalIdx8][interactionVolume.getIndexOnElement(7, 0)].insert(vIdxGlobal);
1700 faceVertices_[globalIdx8][interactionVolume.getIndexOnElement(7, 1)].insert(vIdxGlobal);
1701 faceVertices_[globalIdx8][interactionVolume.getIndexOnElement(7, 2)].insert(vIdxGlobal);
1702 }
1703 }
1704
1705 }
1706}
1707} // end namespace Dumux
1708#endif
Interactionvolume container for 3-d MPFA L-method.
Class including the information of an interaction volume of a MPFA 3D method that does not change wit...
Dune::FieldVector< Scalar, 3 > crossProduct(const Dune::FieldVector< Scalar, 3 > &vec1, const Dune::FieldVector< Scalar, 3 > &vec2)
Cross product of two vectors in three-dimensional Euclidean space.
Definition: math.hh:631
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
Interactionvolume container for 3-d MPFA L-method.
Definition: 3dinteractionvolumecontainer.hh:46
GetPropType< TypeTag, Properties::MPFAInteractionVolume > InteractionVolume
Type for storing an MPFA-interaction-volume. (Usually of type FvMpfaL3dInteractionVolume or FvMpfaL3d...
Definition: 3dinteractionvolumecontainer.hh:94
InteractionVolume & interactionVolume(int vertexIdx)
Returns an interaction volume.
Definition: 3dinteractionvolumecontainer.hh:143
Interactionvolume container for 3-d MPFA L-method on an h-adaptive grid.
Definition: 3dinteractionvolumecontaineradaptive.hh:47
void storeInteractionVolumeInfo()
Stores interaction volumes for each grid vertex.
Definition: 3dinteractionvolumecontaineradaptive.hh:1588
std::set< int > & faceVerticeIndices(int eIdxGlobal, int fIdx)
Returns the set of vertices on an element face.
Definition: 3dinteractionvolumecontaineradaptive.hh:104
FvMpfaL3dInteractionVolumeContainerAdaptive(Problem &problem)
Constructs a FvMpfaL3dInteractionVolumeContainerAdaptive object.
Definition: 3dinteractionvolumecontaineradaptive.hh:113
Properties for a MPFA method.
Base file for properties related to sequential IMPET algorithms.