diff options
Diffstat (limited to 'src/AMRHierLib/AMRGridStitcher.cc')
| -rw-r--r-- | src/AMRHierLib/AMRGridStitcher.cc | 1936 |
1 files changed, 1936 insertions, 0 deletions
diff --git a/src/AMRHierLib/AMRGridStitcher.cc b/src/AMRHierLib/AMRGridStitcher.cc new file mode 100644 index 0000000..61817e4 --- /dev/null +++ b/src/AMRHierLib/AMRGridStitcher.cc @@ -0,0 +1,1936 @@ +//------------------------------------------------------------------------------ +// +// Project: Crack-free Isosurfaces for AMR Data +// Module: $RCSfile$ +// Language: C++ +// Date: $Date$ +// Author: $Author$ +// Version: $Revision$ +// +//------------------------------------------------------------------------------ + +#include "AMRGridStitcher.hh" + +// Own Helper Includes +#include <assume.hh> +#include <Vertex.hh> + +// AMR Specific Includes +#include <AMRGrid.hh> +#include <AMRLevel.hh> + +namespace AMRGridStitcherAuxFct +{ + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism045(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + stitchCellHandler.handleTetrahedronCell(p3, p2, p5, p0); + stitchCellHandler.handlePyramidCell(p1, p4, p5, p2, p0); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism145(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + stitchCellHandler.handleTetrahedronCell(p3, p2, p5, p1); + stitchCellHandler.handlePyramidCell(p0, p3, p5, p4, p1); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism245(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + generateCellsForDeformedTrianglePrism045(stitchCellHandler, p2, p3, p0, p1, p5, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism345(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + generateCellsForDeformedTrianglePrism145(stitchCellHandler, p2, p3, p0, p1, p5, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism0145(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + stitchCellHandler.handleTetrahedronCell(p0, p4, p1, p5); + stitchCellHandler.handlePyramidCell(p0, p1, p2, p3, p5); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism2345(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + generateCellsForDeformedTrianglePrism0145(stitchCellHandler, p2, p3, p0, p1, p5, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism12345(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + stitchCellHandler.handleTetrahedronCell(p0, p4, p1, p3); + stitchCellHandler.handlePyramidCell(p1, p4, p5, p2, p3); + } + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism02345(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + stitchCellHandler.handleTetrahedronCell(p0, p4, p1, p2); + stitchCellHandler.handlePyramidCell(p0, p3, p5, p4, p2); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism01345(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + generateCellsForDeformedTrianglePrism12345(stitchCellHandler, p2, p3, p0, p1, p5, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedTrianglePrism01245(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5) + { + generateCellsForDeformedTrianglePrism02345(stitchCellHandler, p2, p3, p0, p1, p5, p4); + } + template <class StitchCellHandler> + inline void generateCellsForDeformedPyramid04(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4) + { + stitchCellHandler.handleTetrahedronCell(p0, p1, p2, p4); + stitchCellHandler.handleTetrahedronCell(p0, p2, p3, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedPyramid023(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4) + { + stitchCellHandler.handleTetrahedronCell(p0, p1, p2, p4); + stitchCellHandler.handleTetrahedronCell(p0, p2, p3, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedPyramid024(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4) + { + stitchCellHandler.handleTetrahedronCell(p0, p1, p2, p4); + stitchCellHandler.handleTetrahedronCell(p0, p2, p3, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedPyramid1234(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4) + { + stitchCellHandler.handleTetrahedronCell(p0, p1, p3, p4); + stitchCellHandler.handleTetrahedronCell(p1, p2, p3, p4); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedCube024567(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5, const AMRGridCell& p6, const AMRGridCell& p7) + { + generateCellsForDeformedPyramid1234(stitchCellHandler, p3, p7, p4, p0, p2); + generateCellsForDeformedPyramid1234(stitchCellHandler, p1, p0, p4, p5, p2); + stitchCellHandler.handlePyramidCell(p4, p7, p6, p5, p2); + } + + Vertex centroid(const Vertex& v0, const Vertex& v1, const Vertex& v2, const Vertex& v3, const Vertex& v4, const Vertex& v5, const Vertex& v6) + { + return affineCombination(7, &v0, 1.0/7.0, &v1, 1.0/7.0, &v2, 1.0/7.0, &v3, 1.0/7.0, &v4, 1.0/7.0, &v5, 1.0/7.0, &v6, 1.0/7.0); + } + + template <class StitchCellHandler> + inline void generateCellsForDeformedCube1234567(StitchCellHandler& stitchCellHandler, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3, const AMRGridCell& p4, const AMRGridCell& p5, const AMRGridCell& p6, const AMRGridCell& p7) + { + stitchCellHandler.handleTetrahedronCell(p3, p1, p4, p0); + std::cout << "Inserting vertex to avoid illegal tesselation." << std::endl; + Vertex newPos = AMRGridStitcherAuxFct::centroid(p1.pos(), p2.pos(), p3.pos(), p4.pos(), p5.pos(), p6.pos(), p7.pos()); + double newVal = p1.val() + p2.val() + p3.val() + p4.val() + p5.val() + p6.val() + p7.val(); + newVal /= 7; + stitchCellHandler.handleTetrahedronCell(p3, p4, p1, newPos, newVal); + stitchCellHandler.handleTetrahedronCell(p4, p3, p7, newPos, newVal); + stitchCellHandler.handlePyramidCell(p4, p7, p6, p5, newPos, newVal); + stitchCellHandler.handlePyramidCell(p3, p2, p6, p7, newPos, newVal); + stitchCellHandler.handleTetrahedronCell(p1, p2, p3, newPos, newVal); + stitchCellHandler.handlePyramidCell(p1, p5, p6, p2, newPos, newVal); + stitchCellHandler.handleTetrahedronCell(p1, p4, p5, newPos, newVal); + } +} + +template <class StitchCellHandler> +AMRGridStitcher<StitchCellHandler>::AMRGridStitcher(StitchCellHandler& stitchCellHandler, AMRHierarchy& amrHier, AMRGrid *stitchGrid, bool considerNextLevel) +: mStitchCellHandler(stitchCellHandler), mAMRDataSet(amrHier), mStitchGrid(stitchGrid), mConsiderNextLevel(considerNextLevel), mParentLevel(mAMRDataSet.getLevel(stitchGrid->level-1)) +{ + assert(mParentLevel != 0); + // <ATTENTION> Assume spatialrefinement == gridplacementrefinement </ATTENTION> + // <ATTENTION> Assume spatialrefinement % parentLevel->getSpatialRefinemet == 0 </ATTENTION> + int refinementWRTParentLevel[3]; + refinementWRTParentLevel[xAxis] = mStitchGrid->spatialrefinement[xAxis] / mParentLevel->getSpatialRefinement(xAxis); + refinementWRTParentLevel[yAxis] = mStitchGrid->spatialrefinement[yAxis] / mParentLevel->getSpatialRefinement(yAxis); + refinementWRTParentLevel[zAxis] = mStitchGrid->spatialrefinement[zAxis] / mParentLevel->getSpatialRefinement(zAxis); + + mOriginInParentLevel[xAxis] = mStitchGrid->iorigin[xAxis] / refinementWRTParentLevel[xAxis]; + mOriginInParentLevel[yAxis] = mStitchGrid->iorigin[yAxis] / refinementWRTParentLevel[yAxis]; + mOriginInParentLevel[zAxis] = mStitchGrid->iorigin[zAxis] / refinementWRTParentLevel[zAxis]; + // <ATTENTION> Assume refined grid starts at cell boundaries of parent level </ATTENTION> + assume(mStitchGrid->iorigin[xAxis] % refinementWRTParentLevel[xAxis] == 0); + assume(mStitchGrid->iorigin[yAxis] % refinementWRTParentLevel[yAxis] == 0); + assume(mStitchGrid->iorigin[zAxis] % refinementWRTParentLevel[zAxis] == 0); + + mMaxExtInParentLevel[xAxis] = mOriginInParentLevel[xAxis] + (mStitchGrid->dims[xAxis] / refinementWRTParentLevel[xAxis]) -1; + mMaxExtInParentLevel[yAxis] = mOriginInParentLevel[yAxis] + (mStitchGrid->dims[yAxis] / refinementWRTParentLevel[yAxis]) -1; + mMaxExtInParentLevel[zAxis] = mOriginInParentLevel[zAxis] + (mStitchGrid->dims[zAxis] / refinementWRTParentLevel[zAxis]) -1; + // <ATTENTION> Assume refined grid ends at cell boundaries of parent level </ATTENTION> + assume(mStitchGrid->dims[xAxis] % refinementWRTParentLevel[xAxis] == 0); + assume(mStitchGrid->dims[yAxis] % refinementWRTParentLevel[yAxis] == 0); + assume(mStitchGrid->dims[zAxis] % refinementWRTParentLevel[zAxis] == 0); + + // Load all grids necessary + AMRLevel *thisLevel = mAMRDataSet.getLevel(mStitchGrid->level); + for (AMRGrid *currGridPtr = thisLevel->getFirstGrid(); currGridPtr != 0; currGridPtr = currGridPtr->getNextGridInLevel()) + currGridPtr->loadCells(); + + for (AMRGrid *currGridPtr = mParentLevel->getFirstGrid(); currGridPtr != 0; currGridPtr = currGridPtr->getNextGridInLevel()) + currGridPtr->loadCells(); +} + +template <class StitchCellHandler> +AMRGridStitcher<StitchCellHandler>::~AMRGridStitcher() +{ +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::generateStitchCells(bool xyMin, bool xyMax, bool yzMin, bool yzMax, bool xzMin, bool xzMax) +{ + if (xyMin) generateCellsForFace(xAxis, yAxis, zAxis, min); + if (xzMin) generateCellsForFace(zAxis, xAxis, yAxis, min); + if (yzMin) generateCellsForFace(yAxis, zAxis, xAxis, min); + if (xyMax) generateCellsForFace(yAxis, xAxis, zAxis, max); + if (xzMax) generateCellsForFace(xAxis, zAxis, yAxis, max); + if (yzMax) generateCellsForFace(zAxis, yAxis, xAxis, max); + if (xyMin || xzMin) generateCellsForEdge(xAxis, yAxis, min, zAxis, min); + if (xyMax || xzMin) generateCellsForEdge(xAxis, yAxis, min, zAxis, max); + if (xyMin || xzMax) generateCellsForEdge(xAxis, yAxis, max, zAxis, min); + if (xyMax || xzMax) generateCellsForEdge(xAxis, yAxis, max, zAxis, max); + if (xyMin || yzMin) generateCellsForEdge(yAxis, zAxis, min, xAxis, min); + if (xyMin || yzMax) generateCellsForEdge(yAxis, zAxis, min, xAxis, max); + if (xyMax || yzMin) generateCellsForEdge(yAxis, zAxis, max, xAxis, min); + if (xyMax || yzMax) generateCellsForEdge(yAxis, zAxis, max, xAxis, max); + if (xzMin || yzMin) generateCellsForEdge(zAxis, xAxis, min, yAxis, min); + if (xzMin || yzMax) generateCellsForEdge(zAxis, xAxis, min, yAxis, max); + if (xzMax || yzMin) generateCellsForEdge(zAxis, xAxis, max, yAxis, min); + if (xzMax || yzMax) generateCellsForEdge(zAxis, xAxis, max, yAxis, max); + generateCellsForVertex(xAxis, yAxis, zAxis, min, min, min); + generateCellsForVertex(zAxis, yAxis, xAxis, min, min, max); + generateCellsForVertex(xAxis, yAxis, zAxis, max, min, max); + generateCellsForVertex(zAxis, yAxis, xAxis, max, min, min); + generateCellsForVertex(zAxis, yAxis, xAxis, min, max, min); + generateCellsForVertex(xAxis, yAxis, zAxis, max, max, min); + generateCellsForVertex(zAxis, yAxis, xAxis, max, max, max); + generateCellsForVertex(xAxis, yAxis, zAxis, min, max, max); +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectQuadToVertex(const AMRGridCell& qp0, const AMRGridCell& qp1, const AMRGridCell& qp2, const AMRGridCell& qp3, const AMRGridCell& vp, AxisType axis0, AxisType axis1, AxisType axis2, MinMaxType faceType) +{ + // If the coarse grid points are from different grids, i.e. the quad results from several + // coarse grids rather than a single one, only generate stitch cells if we are stitching + // the grid with the smallest index. + if (qp0.containingId() < mStitchGrid->index() || qp1.containingId() < mStitchGrid->index() || + qp2.containingId() < mStitchGrid->index() || qp3.containingId() < mStitchGrid->index()) return; + + if (vp.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= vp.refiningId()) { + // No cell bounadaries are crossed, check only one cell for refinement + AMRGrid *otherGrid = mAMRDataSet.getGridByIndex(vp.refiningId()); + AxisAlignedBox<int> refiningCells = otherGrid->getRefiningCells(vp.grid(), vp.idx(xAxis), vp.idx(yAxis), vp.idx(zAxis)); + int oIdx00[3]; + oIdx00[axis0] = refiningCells.min(axis0); + oIdx00[axis1] = refiningCells.min(axis1); + oIdx00[axis2] = faceType == min ? refiningCells.max(axis2) : refiningCells.min(axis2); + + int oIdx01[3]; + oIdx01[axis0] = oIdx00[axis0] + 1; + oIdx01[axis1] = oIdx00[axis1]; + oIdx01[axis2] = oIdx00[axis2]; + + int oIdx10[3]; + oIdx10[axis0] = oIdx00[axis0]; + oIdx10[axis1] = oIdx00[axis1] + 1; + oIdx10[axis2] = oIdx00[axis2]; + + int oIdx11[3]; + oIdx11[axis0] = oIdx00[axis0] + 1; + oIdx11[axis1] = oIdx00[axis1] + 1; + oIdx11[axis2] = oIdx00[axis2]; + + AMRGridCell vpRef0(otherGrid, oIdx00); + AMRGridCell vpRef1(otherGrid, oIdx10); + AMRGridCell vpRef2(otherGrid, oIdx11); + AMRGridCell vpRef3(otherGrid, oIdx01); + + if (!mConsiderNextLevel || + vpRef0.refiningId() == AMRGrid::NoGrid && + vpRef1.refiningId() == AMRGrid::NoGrid && + vpRef2.refiningId() == AMRGrid::NoGrid && + vpRef3.refiningId() == AMRGrid::NoGrid) + // Only add stitch cell when no stitch cells for next level are generated (i.e. + // one of the refining cells is refined itself by an even finer grid) + mStitchCellHandler.handleHexahedronCell(qp0, qp1, qp2, qp3, vpRef0, vpRef1, vpRef2, vpRef3); + } + // else Do nothing! + } + else { + mStitchCellHandler.handlePyramidCell(qp0, qp1, qp2, qp3, vp); + } +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectQuadToEdge(const AMRGridCell& qp0, const AMRGridCell& qp1, const AMRGridCell& qp2, const AMRGridCell& qp3, const AMRGridCell& ep0, const AMRGridCell& ep1, AxisType axis01, AxisType axis03, AxisType axisPerp0123, MinMaxType qp0Axis01Type, MinMaxType qp0Axis03Type, MinMaxType quadFaceType) +{ + // If the coarse grid points are from different grids, i.e. the quad results from several + // coarse grids rather than a single one, only generate stitch cells if we are stitching + // the grid with the smallest index. + if (qp0.containingId() < mStitchGrid->index() || qp1.containingId() < mStitchGrid->index() || + qp2.containingId() < mStitchGrid->index() || qp3.containingId() < mStitchGrid->index()) return; + + int connectCase = 0; + + AMRGridCell *epRef0 = 0; + AMRGridCell *epRef1 = 0; + if (ep0.refiningId() != AMRGrid::NoGrid) { + connectCase |= 1; + + AMRGrid *otherGrid0 = mAMRDataSet.getGridByIndex(ep0.refiningId()); + AxisAlignedBox<int> refiningCells0 = otherGrid0->getRefiningCells(ep0.grid(), ep0.idx(xAxis), ep0.idx(yAxis), ep0.idx(zAxis)); + + int oIdx0[3]; + oIdx0[axis01] = qp0Axis01Type == min ? refiningCells0.min(axis01) : refiningCells0.max(axis01); + oIdx0[axis03] = qp0Axis03Type == min ? refiningCells0.min(axis03) : refiningCells0.max(axis03); + oIdx0[axisPerp0123] = quadFaceType == min ? refiningCells0.max(axisPerp0123) : refiningCells0.min(axisPerp0123); + epRef0 = new AMRGridCell(otherGrid0, oIdx0); + + int oIdx1[3]; + oIdx1[axis01] = qp0Axis01Type == max ? refiningCells0.min(axis01) : refiningCells0.max(axis01); + oIdx1[axis03] = qp0Axis03Type == min ? refiningCells0.min(axis03) : refiningCells0.max(axis03); + oIdx1[axisPerp0123] = quadFaceType == min ? refiningCells0.max(axisPerp0123) : refiningCells0.min(axisPerp0123); + epRef1 = new AMRGridCell(otherGrid0, oIdx1); + } + + AMRGridCell *epRef2 = 0; + AMRGridCell *epRef3 = 0; + if (ep1.refiningId() != AMRGrid::NoGrid) { + connectCase |= 2; + + AMRGrid *otherGrid1 = mAMRDataSet.getGridByIndex(ep1.refiningId()); + AxisAlignedBox<int> refiningCells1 = otherGrid1->getRefiningCells(ep1.grid(), ep1.idx(xAxis), ep1.idx(yAxis), ep1.idx(zAxis)); + + int oIdx2[3]; + oIdx2[axis01] = qp0Axis01Type == max ? refiningCells1.min(axis01) : refiningCells1.max(axis01); + oIdx2[axis03] = qp0Axis03Type == max ? refiningCells1.min(axis03) : refiningCells1.max(axis03); + oIdx2[axisPerp0123] = quadFaceType == min ? refiningCells1.max(axisPerp0123) : refiningCells1.min(axisPerp0123); + epRef2 = new AMRGridCell(otherGrid1, oIdx2); + + int oIdx3[3]; + oIdx3[axis01] = qp0Axis01Type == min ? refiningCells1.min(axis01) : refiningCells1.max(axis01); + oIdx3[axis03] = qp0Axis03Type == max ? refiningCells1.min(axis03) : refiningCells1.max(axis03); + oIdx3[axisPerp0123] = quadFaceType == min ? refiningCells1.max(axisPerp0123) : refiningCells1.min(axisPerp0123); + epRef3 = new AMRGridCell(otherGrid1, oIdx3); + } + + switch (connectCase) { + case 0: + mStitchCellHandler.handleTrianglePrismCell(qp0, qp1, qp2, qp3, ep0, ep1); + break; + case 1: + assert(epRef0->refiningId() == AMRGrid::NoGrid && + epRef1->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= ep0.refiningId()) { + mStitchCellHandler.handlePyramidCell(*epRef0, *epRef1, qp2, qp3, ep1); + mStitchCellHandler.handleTrianglePrismCell(qp3, qp0, qp1, qp2, *epRef0, *epRef1); + } + break; + case 2: + assert(epRef2->refiningId() == AMRGrid::NoGrid && + epRef3->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= ep1.refiningId()) { + mStitchCellHandler.handlePyramidCell(qp0, qp1, *epRef2, *epRef3, ep0); + mStitchCellHandler.handleTrianglePrismCell(qp1, qp2, qp3, qp0, *epRef2, *epRef3); + } + break; + case 3: + if (!mConsiderNextLevel || + epRef0->refiningId() == AMRGrid::NoGrid && + epRef1->refiningId() == AMRGrid::NoGrid && + epRef2->refiningId() == AMRGrid::NoGrid && + epRef3->refiningId() == AMRGrid::NoGrid) + if (mStitchGrid->index() <= ep0.refiningId() && mStitchGrid->index() <= ep1.refiningId()) { + mStitchCellHandler.handleHexahedronCell(qp0, qp1, qp2, qp3, *epRef0, *epRef1, *epRef2, *epRef3); + } + break; + default: + assert(false); + break; + } + delete epRef0; + delete epRef1; + delete epRef2; + delete epRef3; +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectQuadToQuad(const AMRGridCell& qp0, const AMRGridCell& qp1, const AMRGridCell& qp2, const AMRGridCell& qp3, const AMRGridCell& oqp0, const AMRGridCell& oqp1, const AMRGridCell& oqp2, const AMRGridCell& oqp3, AxisType axis01, AxisType axis03, AxisType axisPerp0123, MinMaxType qp0Axis01Type, MinMaxType qp0Axis03Type, MinMaxType quadFaceType) +{ + // If the coarse grid points are from different grids, i.e. the quad results from several + // coarse grids rather than a single one, only generate stitch cells if we are stitching + // the grid with the smallest index. + if (qp0.containingId() < mStitchGrid->index() || qp1.containingId() < mStitchGrid->index() || + qp2.containingId() < mStitchGrid->index() || qp3.containingId() < mStitchGrid->index()) return; + + + int connectCase = 0; + + AMRGridCell *oqp0Ref = 0; + if (oqp0.refiningId() != AMRGrid::NoGrid) { + connectCase |= 1; + AMRGrid *otherGrid0 = mAMRDataSet.getGridByIndex(oqp0.refiningId()); + AxisAlignedBox<int> refiningCells0 = otherGrid0->getRefiningCells(oqp0.grid(), oqp0.idx(xAxis), oqp0.idx(yAxis), oqp0.idx(zAxis)); + int oIdx0[3]; + oIdx0[axis01] = qp0Axis01Type == max ? refiningCells0.max(axis01) : refiningCells0.min(axis01); + oIdx0[axis03] = qp0Axis03Type == max ? refiningCells0.max(axis03) : refiningCells0.min(axis03); + oIdx0[axisPerp0123] = quadFaceType == min ? refiningCells0.max(axisPerp0123) : refiningCells0.min(axisPerp0123); + oqp0Ref = new AMRGridCell(otherGrid0, oIdx0); + } + + AMRGridCell *oqp1Ref = 0; + if (oqp1.refiningId() != AMRGrid::NoGrid) { + connectCase |= 2; + AMRGrid *otherGrid1 = mAMRDataSet.getGridByIndex(oqp1.refiningId()); + AxisAlignedBox<int> refiningCells1 = otherGrid1->getRefiningCells(oqp1.grid(), oqp1.idx(xAxis), oqp1.idx(yAxis), oqp1.idx(zAxis)); + int oIdx1[3]; + oIdx1[axis01] = qp0Axis01Type == max ? refiningCells1.min(axis01) : refiningCells1.max(axis01); + oIdx1[axis03] = qp0Axis03Type == max ? refiningCells1.max(axis03) : refiningCells1.min(axis03); + oIdx1[axisPerp0123] = quadFaceType == min ? refiningCells1.max(axisPerp0123) : refiningCells1.min(axisPerp0123); + oqp1Ref = new AMRGridCell(otherGrid1, oIdx1); + } + + AMRGridCell *oqp2Ref = 0; + if (oqp2.refiningId() != AMRGrid::NoGrid) { + connectCase |= 4; + AMRGrid *otherGrid2 = mAMRDataSet.getGridByIndex(oqp2.refiningId()); + AxisAlignedBox<int> refiningCells2 = otherGrid2->getRefiningCells(oqp2.grid(), oqp2.idx(xAxis), oqp2.idx(yAxis), oqp2.idx(zAxis)); + int oIdx2[3]; + oIdx2[axis01] = qp0Axis01Type == max ? refiningCells2.min(axis01) : refiningCells2.max(axis01); + oIdx2[axis03] = qp0Axis03Type == max ? refiningCells2.min(axis03) : refiningCells2.max(axis03); + oIdx2[axisPerp0123] = quadFaceType == min ? refiningCells2.max(axisPerp0123) : refiningCells2.min(axisPerp0123); + oqp2Ref = new AMRGridCell(otherGrid2, oIdx2); + } + + AMRGridCell *oqp3Ref = 0; + if (oqp3.refiningId() != AMRGrid::NoGrid) { + connectCase |= 8; + AMRGrid *otherGrid3 = mAMRDataSet.getGridByIndex(oqp3.refiningId()); + AxisAlignedBox<int> refiningCells3 = otherGrid3->getRefiningCells(oqp3.grid(), oqp3.idx(xAxis), oqp3.idx(yAxis), oqp3.idx(zAxis)); + int oIdx3[3]; + oIdx3[axis01] = qp0Axis01Type == max ? refiningCells3.max(axis01) : refiningCells3.min(axis01); + oIdx3[axis03] = qp0Axis03Type == max ? refiningCells3.min(axis03) : refiningCells3.max(axis03); + oIdx3[axisPerp0123] = quadFaceType == min ? refiningCells3.max(axisPerp0123) : refiningCells3.min(axisPerp0123); + oqp3Ref = new AMRGridCell(otherGrid3, oIdx3); + } + + switch (connectCase) { + case 0: + mStitchCellHandler.handleHexahedronCell(qp0, qp1, qp2, qp3, oqp0, oqp1, oqp2, oqp3); + break; + case 1: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp0.refiningId()) { + mStitchCellHandler.handlePyramidCell(qp0, qp1, qp2, qp3, *oqp0Ref); + mStitchCellHandler.handlePyramidCell(qp3, qp2, oqp2, oqp3, *oqp0Ref); + mStitchCellHandler.handlePyramidCell(qp2, qp1, oqp1, oqp2, *oqp0Ref); + } + break; + case 2: + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp1.refiningId()) { + mStitchCellHandler.handlePyramidCell(qp0, qp1, qp2, qp3, *oqp1Ref); + mStitchCellHandler.handlePyramidCell(qp0, qp3, oqp3, oqp0, *oqp1Ref); + mStitchCellHandler.handlePyramidCell(qp3, qp2, oqp2, oqp3, *oqp1Ref); + } + break; + case 3: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp0.refiningId() && mStitchGrid->index() <= oqp1.refiningId()) { + mStitchCellHandler.handleTrianglePrismCell(qp1, qp2, qp3, qp0, *oqp1Ref, *oqp0Ref); + mStitchCellHandler.handleTrianglePrismCell(qp2, oqp2, oqp3, qp3,*oqp1Ref, *oqp0Ref); + } + break; + case 4: + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp2.refiningId()) { + mStitchCellHandler.handlePyramidCell(qp0, qp1, qp2, qp3, *oqp2Ref); + mStitchCellHandler.handlePyramidCell(qp1, qp0, oqp0, oqp1, *oqp2Ref); + mStitchCellHandler.handlePyramidCell(qp0, qp3, oqp3, oqp0, *oqp2Ref); + } + break; + case 5: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp0.refiningId() && mStitchGrid->index() <= oqp2.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube024567(mStitchCellHandler, *oqp0Ref, oqp3, *oqp2Ref, oqp1, qp0, qp3, qp2, qp1); + } + break; + case 6: + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp1.refiningId() && mStitchGrid->index() <= oqp2.refiningId()) { + mStitchCellHandler.handleTrianglePrismCell(qp0, qp1, qp2, qp3, *oqp1Ref, *oqp2Ref); + mStitchCellHandler.handleTrianglePrismCell(oqp0, qp0, qp3, oqp3, *oqp1Ref, *oqp2Ref); + } + break; + case 7: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp2.refiningId() && mStitchGrid->index() <= oqp1.refiningId() && mStitchGrid->index() <= oqp0.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube1234567(mStitchCellHandler, oqp3, *oqp2Ref, *oqp1Ref, *oqp0Ref, qp3, qp2, qp1, qp0); + } + break; + case 8: + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp3.refiningId()) { + mStitchCellHandler.handlePyramidCell(qp0, qp1, qp2, qp3, *oqp3Ref); + mStitchCellHandler.handlePyramidCell(qp1, qp0, oqp0, oqp1, *oqp3Ref); + mStitchCellHandler.handlePyramidCell(qp2, qp1, oqp1, oqp2, *oqp3Ref); + } + break; + case 9: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp0.refiningId() && mStitchGrid->index() <= oqp3.refiningId()) { + mStitchCellHandler.handleTrianglePrismCell(qp0, qp1, qp2, qp3, *oqp0Ref, *oqp3Ref); + mStitchCellHandler.handleTrianglePrismCell(qp1, oqp1, oqp2, qp2, *oqp0Ref, *oqp3Ref); + } + break; + case 10: + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp1.refiningId() && mStitchGrid->index() <= oqp3.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube024567(mStitchCellHandler, *oqp1Ref, oqp0, *oqp3Ref, oqp2, qp1, qp0, qp3, qp2); + } + break; + case 11: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp3.refiningId() && mStitchGrid->index() <= oqp1.refiningId() && mStitchGrid->index() <= oqp0.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube1234567(mStitchCellHandler, oqp2, *oqp1Ref, *oqp0Ref, *oqp3Ref, qp2, qp1, qp0, qp3); + } + break; + case 12: + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp2.refiningId() && mStitchGrid->index() <= oqp3.refiningId()) { + mStitchCellHandler.handleTrianglePrismCell(qp1, qp2, qp3, qp0, *oqp2Ref, *oqp3Ref); + mStitchCellHandler.handleTrianglePrismCell(oqp1, qp1, qp0, oqp0, *oqp2Ref, *oqp3Ref); + } + break; + case 13: + assert(oqp0Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp3.refiningId() && mStitchGrid->index() <= oqp2.refiningId() && mStitchGrid->index() <= oqp0.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube1234567(mStitchCellHandler, oqp1, *oqp0Ref, *oqp3Ref, *oqp2Ref, qp1, qp0, qp3, qp2); + } + break; + case 14: + assert(oqp1Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp2Ref->refiningId() == AMRGrid::NoGrid); + assert(oqp3Ref->refiningId() == AMRGrid::NoGrid); + if (mStitchGrid->index() <= oqp3.refiningId() && mStitchGrid->index() <= oqp2.refiningId() && mStitchGrid->index() <= oqp1.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedCube1234567(mStitchCellHandler, oqp0, *oqp3Ref, *oqp2Ref, *oqp1Ref, qp0, qp3, qp2, qp1); + } + break; + case 15: + if (!mConsiderNextLevel || + oqp0Ref->refiningId() == AMRGrid::NoGrid && + oqp1Ref->refiningId() == AMRGrid::NoGrid && + oqp2Ref->refiningId() == AMRGrid::NoGrid && + oqp3Ref->refiningId() == AMRGrid::NoGrid) + if (mStitchGrid->index() <= oqp0.refiningId() && mStitchGrid->index() <= oqp1.refiningId() && mStitchGrid->index() <= oqp2.refiningId() && mStitchGrid->index() <= oqp3.refiningId()) { + mStitchCellHandler.handleHexahedronCell(qp0, qp1, qp2, qp3, *oqp0Ref, *oqp1Ref, *oqp2Ref, *oqp3Ref); + } + break; + default: + assert(false); + } +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectEdgeSegmentToEdges(const AMRGridCell& ep0, const AMRGridCell& ep1, const AMRGridCell& cp0, const AMRGridCell& cp1, const AMRGridCell& cp2, AxisType edgeAxis, AxisType edgeAxis01, MinMaxType edgeAxis01Type, AxisType edgeAxis02, MinMaxType edgeAxis02Type) +{ + if (ep0.containingId() < mStitchGrid->index() || ep1.containingId() < mStitchGrid->index()) return; + + if (cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + mStitchCellHandler.handleTetrahedronCell(ep0, ep1, cp0, cp1); + mStitchCellHandler.handleTetrahedronCell(ep0, cp0, ep1, cp2); + } + else if (cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0Ref0Idx[3]; + cp0Ref0Idx[edgeAxis] = cp0RefiningCells.min(edgeAxis); + cp0Ref0Idx[edgeAxis01] = edgeAxis01Type == min ? cp0RefiningCells.max(edgeAxis01) : cp0RefiningCells.min(edgeAxis01); + cp0Ref0Idx[edgeAxis02] = edgeAxis02Type == min ? cp0RefiningCells.max(edgeAxis01) : cp0RefiningCells.min(edgeAxis01); + AMRGridCell cp0Ref0(cp0RefiningGrid, cp0Ref0Idx[0], cp0Ref0Idx[1], cp0Ref0Idx[2]); + int cp0Ref1Idx[3]; + cp0Ref1Idx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0Ref1Idx[edgeAxis01] = cp0Ref0Idx[edgeAxis01]; + cp0Ref1Idx[edgeAxis02] = cp0Ref0Idx[edgeAxis02]; + AMRGridCell cp0Ref1(cp0RefiningGrid, cp0Ref1Idx[0], cp0Ref1Idx[1], cp0Ref1Idx[2]); + mStitchCellHandler.handlePyramidCell(ep0, ep1, cp0Ref1, cp0Ref0, cp1); + mStitchCellHandler.handlePyramidCell(ep0, cp0Ref0, cp0Ref1, ep1, cp2); + } + } + else if (cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1Ref0Idx[3]; + cp1Ref0Idx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1Ref0Idx[edgeAxis01] = edgeAxis01Type == min ? cp1RefiningCells.min(edgeAxis01) : cp1RefiningCells.max(edgeAxis01); + cp1Ref0Idx[edgeAxis02] = edgeAxis02Type == min ? cp1RefiningCells.max(edgeAxis02) : cp1RefiningCells.min(edgeAxis02); + AMRGridCell cp1Ref0(cp1RefiningGrid, cp1Ref0Idx[0], cp1Ref0Idx[1], cp1Ref0Idx[2]); + int cp1Ref1Idx[3]; + cp1Ref1Idx[edgeAxis] = cp1RefiningCells.max(edgeAxis); + cp1Ref1Idx[edgeAxis01] = cp1Ref0Idx[edgeAxis01]; + cp1Ref1Idx[edgeAxis02] = cp1Ref0Idx[edgeAxis02]; + AMRGridCell cp1Ref1(cp1RefiningGrid, cp1Ref1Idx[0], cp1Ref1Idx[1], cp1Ref1Idx[2]); + connectQuadToEdge(cp1Ref0, cp1Ref1, ep1, ep0, cp0, cp2, edgeAxis, edgeAxis02, edgeAxis01, min, edgeAxis02Type == min ? max : min, edgeAxis01Type); + } + else { + assert(cp2.refiningId() != AMRGrid::NoGrid); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2Ref0Idx[3]; + cp2Ref0Idx[edgeAxis] = cp2RefiningCells.min(edgeAxis); + cp2Ref0Idx[edgeAxis01] = edgeAxis01Type == min ? cp2RefiningCells.max(edgeAxis01) : cp2RefiningCells.min(edgeAxis01); + cp2Ref0Idx[edgeAxis02] = edgeAxis02Type == min ? cp2RefiningCells.min(edgeAxis02) : cp2RefiningCells.max(edgeAxis02); + AMRGridCell cp2Ref0(cp2RefiningGrid, cp2Ref0Idx[0], cp2Ref0Idx[1], cp2Ref0Idx[2]); + int cp2Ref1Idx[3]; + cp2Ref1Idx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2Ref1Idx[edgeAxis01] = cp2Ref0Idx[edgeAxis01]; + cp2Ref1Idx[edgeAxis02] = cp2Ref0Idx[edgeAxis02]; + AMRGridCell cp2Ref1(cp2RefiningGrid, cp2Ref1Idx[0], cp2Ref1Idx[1], cp2Ref1Idx[2]); + connectQuadToEdge(cp2Ref1, cp2Ref0, ep0, ep1, cp0, cp1, edgeAxis, edgeAxis01, edgeAxis02, max, edgeAxis01Type == min ? max : min, edgeAxis02Type); + } +} + +template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectEdgeSegmentToQuads(const AMRGridCell& ep0, const AMRGridCell& ep1, const AMRGridCell& cp0, const AMRGridCell& cp1, const AMRGridCell& cp2, const AMRGridCell& cp3, const AMRGridCell& cp4, const AMRGridCell& cp5, AxisType edgeAxis, AxisType addAxis0132, MinMaxType edgeAxis0132Type, AxisType addAxis0451, MinMaxType edgeAxis0451Type) +{ + if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 0 + mStitchCellHandler.handleTrianglePrismCell(cp2, cp0, cp1, cp3, ep0, ep1); + mStitchCellHandler.handleTrianglePrismCell(cp0, cp4, cp5, cp1, ep0, ep1); + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 1 + if (mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism145(mStitchCellHandler, cp2, cp0Ref, cp1, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism045(mStitchCellHandler, cp0Ref, cp4, cp5, cp1, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 2 + if (mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism245(mStitchCellHandler, cp2, cp0, cp1Ref, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism345(mStitchCellHandler, cp0, cp4, cp5, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 3 + if (mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + mStitchCellHandler.handleTrianglePrismCell(cp0Ref, ep0, ep1, cp1Ref, cp2, cp3); + mStitchCellHandler.handleTrianglePrismCell(ep0, cp0Ref, cp1Ref, ep1, cp4, cp5); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 4 + if (mStitchGrid->index() <= cp2.refiningId()) { + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + mStitchCellHandler.handleTetrahedronCell(ep1, cp1, cp2Ref, cp3); + std::cout << "Inserting vertex to avoid illegal tesselation (Case 4)." << std::endl; + Vertex newPos = AMRGridStitcherAuxFct::centroid(ep0.pos(), ep1.pos(), cp0.pos(), cp1.pos(), cp2Ref.pos(), cp4.pos(), cp5.pos()); + double newVal = ep0.val() + ep1.val() + cp0.val() + cp1.val() + cp2Ref.val() + cp4.val() + cp5.val(); + newVal /= 7; + mStitchCellHandler.handlePyramidCell(ep0, cp4, cp0, cp2Ref, newPos, newVal); + mStitchCellHandler.handlePyramidCell(ep0, ep1, cp5, cp4, newPos, newVal); + mStitchCellHandler.handlePyramidCell(cp0, cp4, cp5, cp1, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep1, cp1, cp5, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep1, cp2Ref, cp1, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(cp0, cp1, cp2Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp2Ref, ep1, newPos, newVal); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 5 + if (mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism0145(mStitchCellHandler, cp2Ref, cp0Ref, cp1, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism045(mStitchCellHandler, cp0Ref, cp4, cp5, cp1, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 6 + if (mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + mStitchCellHandler.handlePyramidCell(ep0, cp4, cp0, cp2Ref, cp1Ref); + mStitchCellHandler.handlePyramidCell(ep0, ep1, cp5, cp4, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp3, ep1, ep0, cp2Ref, cp1Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 7 + if (mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01245(mStitchCellHandler, cp2Ref, cp0Ref, cp1Ref, cp3, ep0, ep1); + mStitchCellHandler.handleTrianglePrismCell(cp0Ref, cp4, cp5, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 8 + if (mStitchGrid->index() <= cp3.refiningId()) { + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + mStitchCellHandler.handleTetrahedronCell(ep0, cp3Ref, cp0, cp2); + std::cout << "Inserting vertex to avoid illegal tesselation (Case 8)." << std::endl; + Vertex newPos = AMRGridStitcherAuxFct::centroid(ep0.pos(), ep1.pos(), cp0.pos(), cp1.pos(), cp3Ref.pos(), cp4.pos(), cp5.pos()); + double newVal = ep0.val() + ep1.val() + cp0.val() + cp1.val() + cp3Ref.val() + cp4.val() + cp5.val(); + newVal /= 7; + mStitchCellHandler.handlePyramidCell(cp0, cp4, cp5, cp1, newPos, newVal); + mStitchCellHandler.handlePyramidCell(ep0, ep1, cp5, cp4, newPos, newVal); + mStitchCellHandler.handlePyramidCell(ep1, cp3Ref, cp1, cp5, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp4, cp0, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp0, cp3Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(cp0, cp1, cp3Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp3Ref, ep1, newPos, newVal); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 9 + if (mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + + mStitchCellHandler.handlePyramidCell(cp1, cp5, ep1, cp3Ref, cp0Ref); + mStitchCellHandler.handlePyramidCell(ep0, ep1, cp5, cp4, cp0Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp2, cp3Ref, ep1, ep0, cp0Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 10 + if (mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism2345(mStitchCellHandler, cp2, cp0, cp1Ref, cp3Ref, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism345(mStitchCellHandler, cp0, cp4, cp5, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 11 + if (mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism12345(mStitchCellHandler, cp2, cp0Ref, cp1Ref, cp3Ref, ep0, ep1); + mStitchCellHandler.handleTrianglePrismCell(cp0Ref, cp4, cp5, cp1Ref, ep0, ep1); + } + } + else if (!(cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid) && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid) { + // Case 12 ... 15, Case 28 ... 31, Case 44 ... 47 + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + connectQuadToQuad(ep0, cp2Ref, cp3Ref, ep1, cp4, cp0, cp1, cp5, addAxis0451, edgeAxis, addAxis0132, edgeAxis0451Type, max, edgeAxis0132Type); + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 16 + if (mStitchGrid->index() <= cp4.refiningId()) { + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handleTetrahedronCell(ep1, cp4Ref, cp1, cp5); + std::cout << "Inserting vertex to avoid illegal tesselation (Case 16)." << std::endl; + Vertex newPos = AMRGridStitcherAuxFct::centroid(ep0.pos(), ep1.pos(), cp0.pos(), cp1.pos(), cp2.pos(), cp3.pos(), cp4Ref.pos()); + double newVal = ep0.val() + ep1.val() + cp0.val() + cp1.val() + cp2.val() + cp3.val() + cp4Ref.val(); + newVal /= 7; + mStitchCellHandler.handlePyramidCell(ep0, cp4Ref, cp0, cp2, newPos, newVal); + mStitchCellHandler.handlePyramidCell(cp0, cp1, cp3, cp2, newPos, newVal); + mStitchCellHandler.handlePyramidCell(ep0, cp2, cp3, ep1, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, ep1, cp4Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep1, cp1, cp4Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep1, cp3, cp1, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(cp0, cp4Ref, cp1, newPos, newVal); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 17 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism145(mStitchCellHandler, cp2, cp0Ref, cp1, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism0145(mStitchCellHandler, cp0Ref, cp4Ref, cp5, cp1, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 18 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + + mStitchCellHandler.handlePyramidCell(ep0, cp2, cp3, ep1, cp1Ref); + mStitchCellHandler.handlePyramidCell(ep0, cp4Ref, cp0, cp2, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp5, cp4Ref, ep0, ep1, cp1Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 19 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <=cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handleTrianglePrismCell(cp2, cp0Ref, cp1Ref, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01345(mStitchCellHandler, cp0Ref, cp4Ref, cp5, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 20 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp2.refiningId()) { + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp3, ep1, ep0, cp2Ref, cp4Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, cp2Ref, cp0, cp1, cp3, cp4Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, ep1, cp3, cp1, cp5, cp4Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 21 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handlePyramidCell(ep0, cp4Ref, cp0Ref, cp2Ref, ep1); + mStitchCellHandler.handlePyramidCell(cp0Ref, cp4Ref, cp5, cp1, ep1); + mStitchCellHandler.handlePyramidCell(cp0Ref, cp1, cp3, cp2Ref, ep1); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 22 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid024(mStitchCellHandler, ep1, cp3, cp1Ref, cp5, cp4Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid024(mStitchCellHandler, cp1Ref, cp3, cp2Ref, cp0, cp4Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp3, ep1, ep0, cp2Ref, cp4Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 23 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handlePyramidCell(cp0Ref, cp2Ref, ep0, cp4Ref, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp3, ep1, ep0, cp2Ref, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp5, cp4Ref, ep0, ep1, cp1Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 24 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index()) { + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handlePyramidCell(cp0, cp2, ep0, cp4Ref, cp3Ref); + mStitchCellHandler.handleTetrahedronCell(cp0, cp4Ref, cp1, cp3Ref); + mStitchCellHandler.handlePyramidCell(ep1, cp3Ref, cp1, cp5, cp4Ref); + mStitchCellHandler.handleTetrahedronCell(ep0, cp3Ref, ep1, cp4Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 25 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() && cp0.refiningId() <= mStitchGrid->index()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + mStitchCellHandler.handleTrianglePrismCell(cp3Ref, cp1, cp5, ep1, cp0Ref, cp4Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism12345(mStitchCellHandler, cp2, cp3Ref, ep1, ep0, cp0Ref, cp4Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 26 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() && cp1.refiningId() <= mStitchGrid->index()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + + mStitchCellHandler.handleTrianglePrismCell(cp4Ref, cp0, cp2, ep0, cp1Ref, cp3Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism12345(mStitchCellHandler, cp5, cp4Ref, ep0, ep1, cp1Ref, cp3Ref); + } + } + else if (cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 27 + if (mStitchGrid->index() <= cp4.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() && cp1.refiningId() <= mStitchGrid->index() && cp0.refiningId() <= mStitchGrid->index()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism12345(mStitchCellHandler, cp2, cp0Ref, cp1Ref, cp3Ref, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01345(mStitchCellHandler, cp0Ref, cp4Ref, cp5, cp1Ref, ep0, ep1); + } + } + // Case 28 ... 31 see Case 12 ... 15 + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 32 + if (mStitchGrid->index() <= cp5.refiningId()) { + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + mStitchCellHandler.handleTetrahedronCell(ep0, cp0, cp5Ref, cp4); + std::cout << "Inserting vertex to avoid illegal tesselation (Case 32)." << std::endl; + Vertex newPos = AMRGridStitcherAuxFct::centroid(ep0.pos(), ep1.pos(), cp0.pos(), cp1.pos(), cp2.pos(), cp3.pos(), cp5Ref.pos()); + double newVal = ep0.val() + ep1.val() + cp0.val() + cp1.val() + cp2.val() + cp3.val() + cp5Ref.val(); + newVal /= 7; + mStitchCellHandler.handlePyramidCell(ep0, cp2, cp3, ep1, newPos, newVal); + mStitchCellHandler.handlePyramidCell(cp0, cp1, cp3, cp2, newPos, newVal); + mStitchCellHandler.handlePyramidCell(ep1, cp3, cp1, cp5Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, ep1, cp5Ref, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp5Ref, cp0, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(ep0, cp0, cp2, newPos, newVal); + mStitchCellHandler.handleTetrahedronCell(cp0, cp5Ref, cp1, newPos, newVal); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 33 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + mStitchCellHandler.handlePyramidCell(cp1, cp5Ref, ep1, cp3, cp0Ref); + mStitchCellHandler.handlePyramidCell(ep0, cp2, cp3, ep1, cp0Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp4, ep0, ep1, cp5Ref, cp0Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 34 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism245(mStitchCellHandler, cp2, cp0, cp1Ref, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism2345(mStitchCellHandler, cp0, cp4, cp5Ref, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 35 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + mStitchCellHandler.handleTrianglePrismCell(cp2, cp0Ref, cp1Ref, cp3, ep0, ep1); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism02345(mStitchCellHandler, cp0Ref, cp4, cp5Ref, cp1Ref, ep0, ep1); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 36 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp2.refiningId()) { + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism145(mStitchCellHandler, cp3, cp2Ref, cp0, cp1, ep1, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism0145(mStitchCellHandler, cp2Ref, ep0, cp4, cp0, ep1, cp5Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 37 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + mStitchCellHandler.handleTrianglePrismCell(cp3, cp2Ref, cp0Ref, cp1, ep1, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01345(mStitchCellHandler, cp2Ref, ep0, cp4, cp0Ref, ep1, cp5Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 38 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + mStitchCellHandler.handleTrianglePrismCell(ep0, cp4, cp0, cp2Ref, cp5Ref, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01245(mStitchCellHandler, ep1, ep0, cp2Ref, cp3, cp5Ref, cp1Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() != AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 39 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp2.refiningId() && mStitchGrid->index() <= cp1.refiningId() && cp0.refiningId() <= mStitchGrid->index()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp2RefiningGrid = mAMRDataSet.getGridByIndex(cp2.refiningId()); + AxisAlignedBox<int> cp2RefiningCells = cp2RefiningGrid->getRefiningCells(cp2.grid(), cp2.idx(xAxis), cp2.idx(yAxis), cp2.idx(zAxis)); + int cp2RefIdx[3]; + cp2RefIdx[edgeAxis] = cp2RefiningCells.max(edgeAxis); + cp2RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp2RefiningCells.min(addAxis0132) : cp2RefiningCells.max(addAxis0132); + cp2RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp2RefiningCells.max(addAxis0451) : cp2RefiningCells.min(addAxis0451); + AMRGridCell cp2Ref(cp2RefiningGrid, cp2RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism02345(mStitchCellHandler, ep0, cp4, cp0Ref, cp2Ref, cp5Ref, cp1Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedTrianglePrism01245(mStitchCellHandler, ep1, ep0, cp2Ref, cp3, cp5Ref, cp1Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 40 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp3.refiningId()) { + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, ep0, cp4, cp0, cp2, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, cp3Ref, cp2, cp0, cp1, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp2, cp3Ref, ep1, ep0, cp5Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 41 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp3.refiningId() && cp0.refiningId() <= mStitchGrid->index()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid024(mStitchCellHandler, ep0, cp4, cp0Ref, cp2, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid024(mStitchCellHandler, cp3Ref, cp2, cp0Ref, cp1, cp5Ref); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp2, cp3Ref, ep1, ep0, cp5Ref); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + // Case 42 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() <= cp1.refiningId()) { + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + mStitchCellHandler.handlePyramidCell(ep1, cp3Ref, cp1Ref, cp5Ref, ep0); + mStitchCellHandler.handlePyramidCell(cp0, cp4, cp5Ref, cp1Ref, ep0); + mStitchCellHandler.handlePyramidCell(cp0, cp1Ref, cp3Ref, cp2, ep0); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() != AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() != AMRGrid::NoGrid && cp0.refiningId() != AMRGrid::NoGrid) { + // Case 43 + if (mStitchGrid->index() <= cp5.refiningId() && mStitchGrid->index() <= cp3.refiningId() && mStitchGrid->index() <= cp1.refiningId() && mStitchGrid->index() <= cp0.refiningId()) { + AMRGrid *cp0RefiningGrid = mAMRDataSet.getGridByIndex(cp0.refiningId()); + AxisAlignedBox<int> cp0RefiningCells = cp0RefiningGrid->getRefiningCells(cp0.grid(), cp0.idx(xAxis), cp0.idx(yAxis), cp0.idx(zAxis)); + int cp0RefIdx[3]; + cp0RefIdx[edgeAxis] = cp0RefiningCells.max(edgeAxis); + cp0RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp0RefiningCells.max(addAxis0132) : cp0RefiningCells.min(addAxis0132); + cp0RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp0RefiningCells.max(addAxis0451) : cp0RefiningCells.min(addAxis0451); + AMRGridCell cp0Ref(cp0RefiningGrid, cp0RefIdx); + AMRGrid *cp1RefiningGrid = mAMRDataSet.getGridByIndex(cp1.refiningId()); + AxisAlignedBox<int> cp1RefiningCells = cp1RefiningGrid->getRefiningCells(cp1.grid(), cp1.idx(xAxis), cp1.idx(yAxis), cp1.idx(zAxis)); + int cp1RefIdx[3]; + cp1RefIdx[edgeAxis] = cp1RefiningCells.min(edgeAxis); + cp1RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp1RefiningCells.max(addAxis0132) : cp1RefiningCells.min(addAxis0132); + cp1RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp1RefiningCells.max(addAxis0451) : cp1RefiningCells.min(addAxis0451); + AMRGridCell cp1Ref(cp1RefiningGrid, cp1RefIdx); + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[edgeAxis] = cp3RefiningCells.min(edgeAxis); + cp3RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp3RefiningCells.min(addAxis0132) : cp3RefiningCells.max(addAxis0132); + cp3RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp3RefiningCells.max(addAxis0451) : cp3RefiningCells.min(addAxis0451); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + mStitchCellHandler.handlePyramidCell(ep1, cp3Ref, cp1Ref, cp5Ref, ep0); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp4, cp5Ref, cp1Ref, cp0Ref, ep0); + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, cp2, cp0Ref, cp1Ref, cp3Ref, ep0); + } + } + else if (cp5.refiningId() != AMRGrid::NoGrid && cp4.refiningId() != AMRGrid::NoGrid) { + // Case 48 ... 63 + AMRGrid *cp4RefiningGrid = mAMRDataSet.getGridByIndex(cp4.refiningId()); + AxisAlignedBox<int> cp4RefiningCells = cp4RefiningGrid->getRefiningCells(cp4.grid(), cp4.idx(xAxis), cp4.idx(yAxis), cp4.idx(zAxis)); + int cp4RefIdx[3]; + cp4RefIdx[edgeAxis] = cp4RefiningCells.max(edgeAxis); + cp4RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp4RefiningCells.max(addAxis0132) : cp4RefiningCells.min(addAxis0132); + cp4RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp4RefiningCells.min(addAxis0451) : cp4RefiningCells.max(addAxis0451); + AMRGridCell cp4Ref(cp4RefiningGrid, cp4RefIdx); + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[edgeAxis] = cp5RefiningCells.min(edgeAxis); + cp5RefIdx[addAxis0132] = edgeAxis0132Type == min ? cp5RefiningCells.max(addAxis0132) : cp5RefiningCells.min(addAxis0132); + cp5RefIdx[addAxis0451] = edgeAxis0451Type == min ? cp5RefiningCells.min(addAxis0451) : cp5RefiningCells.max(addAxis0451); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + connectQuadToQuad(ep0, ep1, cp5Ref, cp4Ref, cp2, cp3, cp1, cp0, edgeAxis, addAxis0132, addAxis0451, max, edgeAxis0132Type, edgeAxis0451Type); + } + else { + assert(false); + } +} + + template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::connectVertexToQuad(const AMRGridCell& vtx, const AMRGridCell& p0, const AMRGridCell& p1, const AMRGridCell& p2, const AMRGridCell& p3) +{ + assert(p0.refiningId() == AMRGrid::NoGrid); + if (p1.refiningId() == AMRGrid::NoGrid && p2.refiningId() == AMRGrid::NoGrid && p3.refiningId() == AMRGrid::NoGrid) { + mStitchCellHandler.handlePyramidCell(p0, p1, p2, p3, vtx); + } + else if (p1.refiningId() != AMRGrid::NoGrid && p2.refiningId() == AMRGrid::NoGrid && p3.refiningId() == AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p1.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, p1, p2, p3, p0, vtx); + } + } + else if (p1.refiningId() == AMRGrid::NoGrid && p2.refiningId() != AMRGrid::NoGrid && p3.refiningId() == AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p2.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, p2, p3, p0, p1, vtx); + } + } + else if (p1.refiningId() == AMRGrid::NoGrid && p2.refiningId() == AMRGrid::NoGrid && p3.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p3.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid04(mStitchCellHandler, p3, p0, p1, p2, vtx); + } + } + else if (p1.refiningId() != AMRGrid::NoGrid && p2.refiningId() != AMRGrid::NoGrid && p3.refiningId() == AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p1.refiningId() && mStitchGrid->index() <= p2.refiningId()) { + mStitchCellHandler.handlePyramidCell(p0, p1, p2, p3, vtx); + } + } + else if (p1.refiningId() != AMRGrid::NoGrid && p2.refiningId() == AMRGrid::NoGrid && p3.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p1.refiningId() && mStitchGrid->index() <= p3.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid024(mStitchCellHandler, p1, p2, p3, p0, vtx); + } + } + else if (p1.refiningId() == AMRGrid::NoGrid && p2.refiningId() != AMRGrid::NoGrid && p3.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p1.refiningId() && mStitchGrid->index() <= p2.refiningId()) { + mStitchCellHandler.handlePyramidCell(p0, p1, p2, p3, vtx); + } + } + else if (p1.refiningId() != AMRGrid::NoGrid && p2.refiningId() != AMRGrid::NoGrid && p3.refiningId() != AMRGrid::NoGrid) { + if (mStitchGrid->index() <= p1.refiningId() && mStitchGrid->index() <= p2.refiningId() && mStitchGrid->index() <= p3.refiningId()) { + AMRGridStitcherAuxFct::generateCellsForDeformedPyramid1234(mStitchCellHandler, p0, p1, p2, p3, vtx); + } + } + else { + assert(false); + } + +} + + template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::generateCellsForFace(AxisType axis0, AxisType axis1, AxisType axis2, MinMaxType faceType) +{ + assert(axis0 != axis1 && axis0 != axis2 && axis1 != axis2); + + try { + // Here for k=0. k=refiningGrid->dims[2]-1 analogous + // In this step we connect the faces of the refining grid with points in the + // refined grid. + for (int fFaceIdx0=0; fFaceIdx0<mStitchGrid->dims[axis0]-1; ++fFaceIdx0) + for (int fFaceIdx1=0; fFaceIdx1<mStitchGrid->dims[axis1]-1; ++fFaceIdx1) { + int cIdx00[3]; + cIdx00[axis0] = mOriginInParentLevel[axis0] + fFaceIdx0/2; + cIdx00[axis1] = mOriginInParentLevel[axis1] + fFaceIdx1/2; + if (faceType == min) + cIdx00[axis2] = mOriginInParentLevel[axis2] - 1; + else + cIdx00[axis2] = mMaxExtInParentLevel[axis2] + 1; + AMRGridCell cp0 = mParentLevel->getGridCell(cIdx00); + + int cIdx01[3]; + cIdx01[axis0] = cIdx00[axis0] + 1; + cIdx01[axis1] = cIdx00[axis1]; + cIdx01[axis2] = cIdx00[axis2]; + AMRGridCell cp3 = mParentLevel->getGridCell(cIdx01); + + int cIdx10[3]; + cIdx10[axis0] = cIdx00[axis0]; + cIdx10[axis1] = cIdx00[axis1] + 1; + cIdx10[axis2] = cIdx00[axis2]; + AMRGridCell cp1 = mParentLevel->getGridCell(cIdx10); + + int cIdx11[3]; + cIdx11[axis0] = cIdx00[axis0] + 1; + cIdx11[axis1] = cIdx00[axis1] + 1; + cIdx11[axis2] = cIdx00[axis2]; + AMRGridCell cp2 = mParentLevel->getGridCell(cIdx11);; + + int fIdx00[3]; + fIdx00[axis0] = fFaceIdx0; + fIdx00[axis1] = fFaceIdx1; + if (faceType == min) + fIdx00[axis2] = 0; + else + fIdx00[axis2] = mStitchGrid->dims[axis2] - 1; + AMRGridCell fp0(mStitchGrid, fIdx00[0], fIdx00[1], fIdx00[2]); + + int fIdx01[3]; + fIdx01[axis0] = fIdx00[axis0] + 1; + fIdx01[axis1] = fIdx00[axis1]; + fIdx01[axis2] = fIdx00[axis2]; + AMRGridCell fp3(mStitchGrid, fIdx01[0], fIdx01[1], fIdx01[2]); + + int fIdx10[3]; + fIdx10[axis0] = fIdx00[axis0]; + fIdx10[axis1] = fIdx00[axis1] + 1; + fIdx10[axis2] = fIdx00[axis2]; + AMRGridCell fp1(mStitchGrid, fIdx10[0], fIdx10[1], fIdx10[2]); + + int fIdx11[3]; + fIdx11[axis0] = fIdx00[axis0] + 1; + fIdx11[axis1] = fIdx00[axis1] + 1; + fIdx11[axis2] = fIdx00[axis2]; + AMRGridCell fp2(mStitchGrid, fIdx11[0], fIdx11[1], fIdx11[2]); + + // Compute combination from 2D case. If i is even, we connect with the corresponding + // point, if it is odd, we connect with the corresponding edge. The same for j. + // The combination yields: + // Axis 1 Axis 0 Result + // ====================== + // point point point + // point edge edge 0 (parallel axis 0) + // edge point edge 1 (parallel to axis 1) + // edge edge quad + enum ConnectionTypeType { Point = 0, Edge0 = 1, Edge1 = 2, Quad = 3 }; + ConnectionTypeType connectionType = ConnectionTypeType((fFaceIdx1%2) << 1 | fFaceIdx0%2); + switch (connectionType) { + case Point: + connectQuadToVertex(fp0, fp1, fp2, fp3, cp0, axis0, axis1, axis2, faceType); + break; + case Edge0: + //connectQuadToAxis0Edge(fp0, fp1, fp2, fp3, cp0, cp3, axis0, axis1, axis2, faceType); + connectQuadToEdge(fp0, fp1, fp2, fp3, cp0, cp3, axis1, axis0, axis2, min, max, faceType); + break; + case Edge1: + //connectQuadToAxis1Edge(fp0, fp1, fp2, fp3, cp0, cp1, axis0, axis1, axis2, faceType); + connectQuadToEdge(fp3, fp0, fp1, fp2, cp0, cp1, axis0, axis1, axis2, max, max, faceType); + break; + case Quad: + connectQuadToQuad(fp0, fp1, fp2, fp3, cp0, cp1, cp2, cp3, axis1, axis0, axis2, max, max, faceType); + break; + default: + std::cerr << "Internal Error (generateCellsForFace(AxisType, AxisType, AxisType)): Invalid enum value." << std::endl; + abort(); + } + } + } + catch (AMRLevel::NoGrid ng) { + std::cerr << "Error (AMRGridStitcher::generateCellsForFace(AxisType, AxisType, AxisType, MinMaxType): " << ng << std::endl; + } +} + + template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::generateCellsForEdge(AxisType edgeAxis, AxisType nonEdgeAxis0, MinMaxType nonEdgeAxis0Type, AxisType nonEdgeAxis1, MinMaxType nonEdgeAxis1Type) +{ + AxisType addAxis0132; + AxisType addAxis0451; + MinMaxType addAxis0132Type; + MinMaxType addAxis0451Type; + if (nonEdgeAxis1Type == min && nonEdgeAxis0Type == min) { + addAxis0132 = nonEdgeAxis0; + addAxis0451 = nonEdgeAxis1; + addAxis0132Type = min; + addAxis0451Type = min; + } + else if (nonEdgeAxis1Type == min && nonEdgeAxis0Type == max) { + addAxis0132 = nonEdgeAxis1; + addAxis0451 = nonEdgeAxis0; + addAxis0132Type = min; + addAxis0451Type = max; + } + if (nonEdgeAxis1Type == max && nonEdgeAxis0Type == min) { + addAxis0132 = nonEdgeAxis1; + addAxis0451 = nonEdgeAxis0; + addAxis0132Type = max; + addAxis0451Type = min; + } + if (nonEdgeAxis1Type == max && nonEdgeAxis0Type == max) { + addAxis0132 = nonEdgeAxis0; + addAxis0451 = nonEdgeAxis1; + addAxis0132Type = max; + addAxis0451Type = max; + } + + try { + for (int edgeIdx=0; edgeIdx<mStitchGrid->dims[edgeAxis]-1; ++edgeIdx) { + int cIdx0[3]; + cIdx0[edgeAxis] = mOriginInParentLevel[edgeAxis] + edgeIdx/2; + if (nonEdgeAxis0Type == min) + cIdx0[nonEdgeAxis0] = mOriginInParentLevel[nonEdgeAxis0] - 1; + else + cIdx0[nonEdgeAxis0] = mMaxExtInParentLevel[nonEdgeAxis0] + 1; + if (nonEdgeAxis1Type == min) + cIdx0[nonEdgeAxis1] = mOriginInParentLevel[nonEdgeAxis1] - 1; + else + cIdx0[nonEdgeAxis1] = mMaxExtInParentLevel[nonEdgeAxis1] + 1; + AMRGridCell cp0 = mParentLevel->getGridCell(cIdx0); + + int cIdx1[3]; + cIdx1[edgeAxis] = cIdx0[edgeAxis] + 1; + cIdx1[nonEdgeAxis0] = cIdx0[nonEdgeAxis0]; + cIdx1[nonEdgeAxis1] = cIdx0[nonEdgeAxis1]; + AMRGridCell cp1 = mParentLevel->getGridCell(cIdx1); + + int cIdx2[3]; + cIdx2[edgeAxis] = cIdx0[edgeAxis]; + if (addAxis0132Type == min) + cIdx2[addAxis0132] = cIdx0[addAxis0132] + 1; + else + cIdx2[addAxis0132] = cIdx0[addAxis0132] - 1; + cIdx2[addAxis0451] = cIdx0[addAxis0451]; + AMRGridCell cp2 = mParentLevel->getGridCell(cIdx2); + + int cIdx3[3]; + cIdx3[edgeAxis] = cIdx2[edgeAxis] + 1; + cIdx3[addAxis0132] = cIdx2[addAxis0132]; + cIdx3[addAxis0451] = cIdx2[addAxis0451]; + AMRGridCell cp3 = mParentLevel->getGridCell(cIdx3); + + int cIdx4[3]; + cIdx4[edgeAxis] = cIdx0[edgeAxis]; + cIdx4[addAxis0132] = cIdx0[addAxis0132]; + if (addAxis0451Type == min) + cIdx4[addAxis0451] = cIdx0[addAxis0451] + 1; + else + cIdx4[addAxis0451] = cIdx0[addAxis0451] - 1; + AMRGridCell cp4 = mParentLevel->getGridCell(cIdx4); + + int cIdx5[3]; + cIdx5[edgeAxis] = cIdx4[edgeAxis] + 1; + cIdx5[addAxis0132] = cIdx4[addAxis0132]; + cIdx5[addAxis0451] = cIdx4[addAxis0451]; + AMRGridCell cp5 = mParentLevel->getGridCell(cIdx5); + + int eIdx0[3]; + eIdx0[edgeAxis] = edgeIdx; + if (nonEdgeAxis0Type == min) + eIdx0[nonEdgeAxis0] = 0; + else + eIdx0[nonEdgeAxis0] = mStitchGrid->dims[nonEdgeAxis0]-1; + if (nonEdgeAxis1Type == min) + eIdx0[nonEdgeAxis1] = 0; + else + eIdx0[nonEdgeAxis1] = mStitchGrid->dims[nonEdgeAxis1]-1; + AMRGridCell ep0(mStitchGrid, eIdx0[0], eIdx0[1], eIdx0[2]); + + int eIdx1[3]; + eIdx1[edgeAxis] = eIdx0[edgeAxis] + 1; + eIdx1[nonEdgeAxis0] = eIdx0[nonEdgeAxis0]; + eIdx1[nonEdgeAxis1] = eIdx0[nonEdgeAxis1]; + AMRGridCell ep1(mStitchGrid, eIdx1[0], eIdx1[1], eIdx1[2]); + + if (edgeIdx % 2) { + connectEdgeSegmentToQuads(ep0, ep1, cp0, cp1, cp2, cp3, cp4, cp5, edgeAxis, addAxis0132, addAxis0132Type, addAxis0451, addAxis0451Type); + } + else { + connectEdgeSegmentToEdges(ep0, ep1, cp0, cp2, cp4, edgeAxis, addAxis0132, addAxis0132Type, addAxis0451, addAxis0451Type); + } + } + } + catch (AMRLevel::NoGrid ng) { + std::cerr << "Error (AMRGridStitcher::generateCellsForEdge(AxisType, AxisType, MinMaxType, AxisType, MinMaxType): " << ng << std::endl; + } +} + + template<class StitchCellHandler> +void AMRGridStitcher<StitchCellHandler>::generateCellsForVertex(AxisType axis0, AxisType axis1, AxisType axis2, MinMaxType axis0VertexType, MinMaxType axis1VertexType, MinMaxType axis2VertexType) +{ + try { + int fIdx[3]; + fIdx[axis0] = axis0VertexType == min ? 0 : mStitchGrid->dims[axis0] - 1; + fIdx[axis1] = axis1VertexType == min ? 0 : mStitchGrid->dims[axis1] - 1; + fIdx[axis2] = axis2VertexType == min ? 0 : mStitchGrid->dims[axis2] - 1; + AMRGridCell fp(mStitchGrid, fIdx); + + int cp0Idx[3]; + cp0Idx[axis0] = axis0VertexType == min ? mOriginInParentLevel[axis0] - 1 : mMaxExtInParentLevel[axis0] + 1; + cp0Idx[axis1] = axis1VertexType == min ? mOriginInParentLevel[axis1] - 1 : mMaxExtInParentLevel[axis1] + 1; + cp0Idx[axis2] = axis2VertexType == min ? mOriginInParentLevel[axis2] - 1 : mMaxExtInParentLevel[axis2] + 1; + AMRGridCell cp0 = mParentLevel->getGridCell(cp0Idx); + + int cp1Idx[3]; + cp1Idx[axis0] = cp0Idx[axis0]; + cp1Idx[axis1] = axis1VertexType == min ? cp0Idx[axis1] + 1 : cp0Idx[axis1] - 1; + cp1Idx[axis2] = cp0Idx[axis2]; + AMRGridCell cp1 = mParentLevel->getGridCell(cp1Idx); + + int cp2Idx[3]; + cp2Idx[axis0] = axis0VertexType == min ? cp0Idx[axis0] + 1 : cp0Idx[axis0] - 1; + cp2Idx[axis1] = cp0Idx[axis1]; + cp2Idx[axis2] = cp0Idx[axis2]; + AMRGridCell cp2 = mParentLevel->getGridCell(cp2Idx); + + int cp3Idx[3]; + cp3Idx[axis0] = axis0VertexType == min ? cp1Idx[axis0] + 1 : cp1Idx[axis0] - 1; + cp3Idx[axis1] = cp1Idx[axis1]; + cp3Idx[axis2] = cp1Idx[axis2]; + AMRGridCell cp3 = mParentLevel->getGridCell(cp3Idx); + + int cp4Idx[3]; + cp4Idx[axis0] = cp0Idx[axis0]; + cp4Idx[axis1] = cp0Idx[axis1]; + cp4Idx[axis2] = axis2VertexType == min ? cp0Idx[axis2] + 1 : cp0Idx[axis2] - 1; + AMRGridCell cp4 = mParentLevel->getGridCell(cp4Idx); + + int cp5Idx[3]; + cp5Idx[axis0] = cp1Idx[axis0]; + cp5Idx[axis1] = cp1Idx[axis1]; + cp5Idx[axis2] = axis2VertexType == min ? cp1Idx[axis2] + 1 : cp1Idx[axis2] - 1; + AMRGridCell cp5 = mParentLevel->getGridCell(cp5Idx); + + int cp6Idx[3]; + cp6Idx[axis0] = cp2Idx[axis0]; + cp6Idx[axis1] = cp0Idx[axis1]; + cp6Idx[axis2] = cp4Idx[axis2]; + AMRGridCell cp6 = mParentLevel->getGridCell(cp6Idx); + + if (cp6.refiningId() == AMRGrid::NoGrid && cp5.refiningId() == AMRGrid::NoGrid && cp4.refiningId() == AMRGrid::NoGrid && cp3.refiningId() == AMRGrid::NoGrid && cp2.refiningId() == AMRGrid::NoGrid && cp1.refiningId() == AMRGrid::NoGrid && cp0.refiningId() == AMRGrid::NoGrid) { + mStitchCellHandler.handlePyramidCell(cp0, cp2, cp3, cp1, fp); + mStitchCellHandler.handlePyramidCell(cp0, cp1, cp5, cp4, fp); + mStitchCellHandler.handlePyramidCell(cp0, cp4, cp6, cp2, fp); + } + else if (!(cp3.refiningId() != AMRGrid::NoGrid || cp5.refiningId() != AMRGrid::NoGrid || cp6.refiningId() != AMRGrid::NoGrid)) { + connectVertexToQuad(fp, cp3, cp1, cp0, cp2); + connectVertexToQuad(fp, cp5, cp4, cp0, cp1); + connectVertexToQuad(fp, cp6, cp2, cp0, cp4); + } + else if (cp3.refiningId() != AMRGrid::NoGrid) { + AMRGrid *cp3RefiningGrid = mAMRDataSet.getGridByIndex(cp3.refiningId()); + AxisAlignedBox<int> cp3RefiningCells = cp3RefiningGrid->getRefiningCells(cp3.grid(), cp3.idx(xAxis), cp3.idx(yAxis), cp3.idx(zAxis)); + int cp3RefIdx[3]; + cp3RefIdx[axis0] = axis0VertexType == min ? cp3RefiningCells.min(axis0) : cp3RefiningCells.max(axis0); + cp3RefIdx[axis1] = axis1VertexType == min ? cp3RefiningCells.min(axis1) : cp3RefiningCells.max(axis1); + cp3RefIdx[axis2] = axis2VertexType == min ? cp3RefiningCells.max(axis2) : cp3RefiningCells.min(axis2); + AMRGridCell cp3Ref(cp3RefiningGrid, cp3RefIdx); + if (axis2VertexType == min) + connectEdgeSegmentToQuads(cp3Ref, fp, cp0, cp4, cp2, cp6, cp1, cp5, axis2, axis0, axis0VertexType, axis1, axis1VertexType); + else + connectEdgeSegmentToQuads(fp, cp3Ref, cp4, cp0, cp5, cp1, cp6, cp2, axis2, axis1, axis1VertexType, axis0, axis0VertexType); + } + else if (cp5.refiningId() != AMRGrid::NoGrid) { + AMRGrid *cp5RefiningGrid = mAMRDataSet.getGridByIndex(cp5.refiningId()); + AxisAlignedBox<int> cp5RefiningCells = cp5RefiningGrid->getRefiningCells(cp5.grid(), cp5.idx(xAxis), cp5.idx(yAxis), cp5.idx(zAxis)); + int cp5RefIdx[3]; + cp5RefIdx[axis0] = axis0VertexType == min ? cp5RefiningCells.max(axis0) : cp5RefiningCells.min(axis0); + cp5RefIdx[axis1] = axis1VertexType == min ? cp5RefiningCells.min(axis1) : cp5RefiningCells.max(axis1); + cp5RefIdx[axis2] = axis2VertexType == min ? cp5RefiningCells.min(axis2) : cp5RefiningCells.max(axis2); + AMRGridCell cp5Ref(cp5RefiningGrid, cp5RefIdx); + if (axis0VertexType == min) + connectEdgeSegmentToQuads(cp5Ref, fp, cp0, cp2, cp1, cp3, cp4, cp6, axis0, axis1, axis1VertexType, axis2, axis2VertexType); + else + connectEdgeSegmentToQuads(fp, cp5Ref, cp2, cp0, cp6, cp4, cp3, cp1, axis0, axis2, axis2VertexType, axis1, axis1VertexType); + } + else if (cp6.refiningId() != AMRGrid::NoGrid) { + AMRGrid *cp6RefiningGrid = mAMRDataSet.getGridByIndex(cp6.refiningId()); + AxisAlignedBox<int> cp6RefiningCells = cp6RefiningGrid->getRefiningCells(cp6.grid(), cp6.idx(xAxis), cp6.idx(yAxis), cp6.idx(zAxis)); + int cp6RefIdx[3]; + cp6RefIdx[axis0] = axis0VertexType == min ? cp6RefiningCells.min(axis0) : cp6RefiningCells.max(axis0); + cp6RefIdx[axis1] = axis1VertexType == min ? cp6RefiningCells.max(axis1) : cp6RefiningCells.min(axis1); + cp6RefIdx[axis2] = axis2VertexType == min ? cp6RefiningCells.min(axis2) : cp6RefiningCells.max(axis2); + AMRGridCell cp6Ref(cp6RefiningGrid, cp6RefIdx); + if (axis0VertexType == min) + connectEdgeSegmentToQuads(cp6Ref, fp, cp0, cp1, cp4, cp5, cp2, cp3, axis1, axis2, axis2VertexType, axis0, axis0VertexType); + else + connectEdgeSegmentToQuads(fp, cp6Ref, cp1, cp0, cp3, cp2, cp5, cp4, axis1, axis0, axis0VertexType, axis2, axis2VertexType); + } + } + catch (AMRLevel::NoGrid ng) { + std::cerr << "Error (AMRGridStitcher::generateCellsForVertex(AxisType, AxisType, AxisType, MinMaxType, MinMaxType, MinMaxType): " << ng << std::endl; + } +} |