diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H index 6dd0675993..9a22055dc9 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H @@ -600,6 +600,19 @@ private: PackedBoolList& boundaryFacesToRemove ); + void calcNeighbourCellCentres + ( + const polyMesh& mesh, + const pointField& cellCentres, + pointField& neiCc + ) const; + + void selectSeparatedCoupledFaces + ( + const polyMesh& mesh, + boolList& selected + ) const; + //- From meshRefinementBaffles.C. Use insidePoint for a surface to // determine the cell zone. void findCellZoneInsideWalk diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C index c943ca978b..4e70aab964 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C @@ -2744,12 +2744,45 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches } else { - // internal face - faces[dualFaceI] = newDualFace; - owner[dualFaceI] = own; - neighbour[dualFaceI] = nei; +// if +// ( +// ptPairs_.isPointPair(vA, vB) +// || ftPtConformer_.featurePointPairs().isPointPair(vA, vB) +// ) +// { + patchIndex = geometryToConformTo_.findPatch(ptA, ptB); +// } - dualFaceI++; + if (patchIndex != -1) + { +// patchFaces[patchIndex].append(newDualFace); +// patchOwners[patchIndex].append(own); +// indirectPatchFace[patchIndex].append(false); +// +// if +// ( +// labelPair(vB->index(), vB->procIndex()) +// < labelPair(vA->index(), vA->procIndex()) +// ) +// { +// patchPPSlaves[patchIndex].append(vB->index()); +// } +// else +// { +// patchPPSlaves[patchIndex].append(vA->index()); +// } + +// baffleFaces[dualFaceI] = patchIndex; + } +// else + { + // internal face + faces[dualFaceI] = newDualFace; + owner[dualFaceI] = own; + neighbour[dualFaceI] = nei; + + dualFaceI++; + } } } } diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C index 07fd509e1f..0de3572e00 100644 --- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C +++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C @@ -111,7 +111,7 @@ void Foam::conformalVoronoiMesh::writeMesh(const fileName& instance) { pointField points; - labelList boundaryPts(number_of_finite_cells(), -1); + labelList boundaryPts; faceList faces; labelList owner; labelList neighbour; @@ -407,6 +407,78 @@ void Foam::conformalVoronoiMesh::writeMesh(const fileName& instance) } +void Foam::conformalVoronoiMesh::calcNeighbourCellCentres +( + const polyMesh& mesh, + const pointField& cellCentres, + pointField& neiCc +) const +{ + label nBoundaryFaces = mesh.nFaces() - mesh.nInternalFaces(); + + if (neiCc.size() != nBoundaryFaces) + { + FatalErrorIn("conformalVoronoiMesh::calcNeighbourCellCentres(..)") + << "nBoundaries:" << nBoundaryFaces + << " neiCc:" << neiCc.size() + << abort(FatalError); + } + + const polyBoundaryMesh& patches = mesh.boundaryMesh(); + + forAll(patches, patchI) + { + const polyPatch& pp = patches[patchI]; + + const labelUList& faceCells = pp.faceCells(); + + label bFaceI = pp.start() - mesh.nInternalFaces(); + + if (pp.coupled()) + { + forAll(faceCells, i) + { + neiCc[bFaceI] = cellCentres[faceCells[i]]; + bFaceI++; + } + } + } + + // Swap coupled boundaries. Apply separation to cc since is coordinate. + syncTools::swapBoundaryFacePositions(mesh, neiCc); +} + + +void Foam::conformalVoronoiMesh::selectSeparatedCoupledFaces +( + const polyMesh& mesh, + boolList& selected +) const +{ + const polyBoundaryMesh& patches = mesh.boundaryMesh(); + + forAll(patches, patchI) + { + // Check all coupled. Avoid using .coupled() so we also pick up AMI. + if (isA(patches[patchI])) + { + const coupledPolyPatch& cpp = refCast + ( + patches[patchI] + ); + + if (cpp.separated() || !cpp.parallel()) + { + forAll(cpp, i) + { + selected[cpp.start()+i] = true; + } + } + } + } +} + + void Foam::conformalVoronoiMesh::findCellZoneInsideWalk ( const polyMesh& mesh, @@ -417,7 +489,7 @@ void Foam::conformalVoronoiMesh::findCellZoneInsideWalk { // Analyse regions. Reuse regionsplit boolList blockedFace(mesh.nFaces()); - //selectSeparatedCoupledFaces(blockedFace); + selectSeparatedCoupledFaces(mesh, blockedFace); forAll(faceToSurface, faceI) { @@ -630,42 +702,90 @@ void Foam::conformalVoronoiMesh::calcFaceZones const labelList& faceOwner = mesh.faceOwner(); const labelList& faceNeighbour = mesh.faceNeighbour(); + labelList neiFaceOwner(mesh.nFaces() - mesh.nInternalFaces(), -1); + + const polyBoundaryMesh& patches = mesh.boundaryMesh(); + + forAll(patches, patchI) + { + const polyPatch& pp = patches[patchI]; + + const labelUList& faceCells = pp.faceCells(); + + label bFaceI = pp.start() - mesh.nInternalFaces(); + + if (pp.coupled()) + { + forAll(faceCells, i) + { + neiFaceOwner[bFaceI] = cellToSurface[faceCells[i]]; + bFaceI++; + } + } + } + + syncTools::swapBoundaryFaceList(mesh, neiFaceOwner); + forAll(faces, faceI) { const label ownerSurfaceI = cellToSurface[faceOwner[faceI]]; + if (faceToSurface[faceI] >= 0) + { + continue; + } + if (mesh.isInternalFace(faceI)) { const label neiSurfaceI = cellToSurface[faceNeighbour[faceI]]; - flipMap[faceI] = - ( - ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI) - ? false - : true - ); - if ( (ownerSurfaceI >= 0 || neiSurfaceI >= 0) && ownerSurfaceI != neiSurfaceI ) { - if (ownerSurfaceI > neiSurfaceI) - { - faceToSurface[faceI] = ownerSurfaceI; - } - else - { - faceToSurface[faceI] = neiSurfaceI; - } + flipMap[faceI] = + ( + ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI) + ? false + : true + ); + + faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI); } } else { - if (ownerSurfaceI >= 0) + label patchID = mesh.boundaryMesh().whichPatch(faceI); + + if (mesh.boundaryMesh()[patchID].coupled()) { - faceToSurface[faceI] = ownerSurfaceI; + const label neiSurfaceI = + neiFaceOwner[faceI - mesh.nInternalFaces()]; + + if + ( + (ownerSurfaceI >= 0 || neiSurfaceI >= 0) + && ownerSurfaceI != neiSurfaceI + ) + { + flipMap[faceI] = + ( + ownerSurfaceI == max(ownerSurfaceI, neiSurfaceI) + ? false + : true + ); + + faceToSurface[faceI] = max(ownerSurfaceI, neiSurfaceI); + } + } + else + { + if (ownerSurfaceI >= 0) + { + faceToSurface[faceI] = ownerSurfaceI; + } } } } @@ -674,104 +794,150 @@ void Foam::conformalVoronoiMesh::calcFaceZones const PtrList& surfZones = geometryToConformTo().surfZones(); - labelList insidePointNamedSurfaces + labelList unclosedSurfaces ( - surfaceZonesInfo::getInsidePointNamedSurfaces(surfZones) + surfaceZonesInfo::getUnclosedNamedSurfaces + ( + surfZones, + geometryToConformTo().geometry(), + geometryToConformTo().surfaces() + ) ); + pointField neiCc(mesh.nFaces() - mesh.nInternalFaces()); + calcNeighbourCellCentres + ( + mesh, + cellCentres, + neiCc + ); + + OBJstream intersections(time().path()/"ints.obj"); + OBJstream intersectionFaces(time().path()/"intFaces.obj"); + // Use intersection of cellCentre connections forAll(faces, faceI) { - if - ( - mesh.isInternalFace(faceI) - && faceToSurface[faceI] < 0 - ) + if (faceToSurface[faceI] >= 0) + { + continue; + } + + label patchID = mesh.boundaryMesh().whichPatch(faceI); + + const label own = faceOwner[faceI]; + + List surfHit; + labelList hitSurface; + + if (mesh.isInternalFace(faceI)) { - const label own = faceOwner[faceI]; const label nei = faceNeighbour[faceI]; - pointIndexHit surfHit; - label hitSurface; - - geometryToConformTo().findSurfaceAnyIntersection + geometryToConformTo().findSurfaceAllIntersections ( cellCentres[own], cellCentres[nei], surfHit, hitSurface ); + } + else if (patchID != -1 && mesh.boundaryMesh()[patchID].coupled()) + { + geometryToConformTo().findSurfaceAllIntersections + ( + cellCentres[own], + neiCc[faceI - mesh.nInternalFaces()], + surfHit, + hitSurface + ); - if (surfHit.hit()) + if (surfHit.size() == 1 && surfHit[0].hit()) { - if (findIndex(insidePointNamedSurfaces, hitSurface) != -1) - { - faceToSurface[faceI] = hitSurface; - - vectorField norm; - geometryToConformTo().getNormal + intersections.write + ( + linePointRef ( - hitSurface, - List(1, surfHit), - norm - ); - - const vector fN = faces[faceI].normal(mesh.points()); - - if ((norm[0] & fN/(mag(fN) + SMALL)) < 0) - { - flipMap[faceI] = true; - } - else - { - flipMap[faceI] = false; - } - } + cellCentres[own], + neiCc[faceI - mesh.nInternalFaces()] + ) + ); } } - } - - labelList neiCellSurface(mesh.nFaces()-mesh.nInternalFaces()); - const polyBoundaryMesh& patches = mesh.boundaryMesh(); - - forAll(patches, patchI) - { - const polyPatch& pp = patches[patchI]; - - if (pp.coupled()) + // If there are multiple intersections then do not add to + // a faceZone + if (surfHit.size() == 1 && surfHit[0].hit()) { - forAll(pp, i) + if (findIndex(unclosedSurfaces, hitSurface[0]) != -1) { - label faceI = pp.start()+i; - label ownSurface = cellToSurface[faceOwner[faceI]]; - neiCellSurface[faceI - mesh.nInternalFaces()] = ownSurface; - } - } - } - syncTools::swapBoundaryFaceList(mesh, neiCellSurface); + vectorField norm; + geometryToConformTo().getNormal + ( + hitSurface[0], + List(1, surfHit[0]), + norm + ); - forAll(patches, patchI) - { - const polyPatch& pp = patches[patchI]; + vector fN = faces[faceI].normal(mesh.points()); + fN /= mag(fN) + SMALL; - if (pp.coupled()) - { - forAll(pp, i) - { - label faceI = pp.start()+i; - label ownSurface = cellToSurface[faceOwner[faceI]]; - label neiSurface = neiCellSurface[faceI-mesh.nInternalFaces()]; - - if (faceToSurface[faceI] == -1 && (ownSurface != neiSurface)) + if ((norm[0] & fN) < 0) { - // Give face the max cell zone - faceToSurface[faceI] = max(ownSurface, neiSurface); + flipMap[faceI] = true; } + else + { + flipMap[faceI] = false; + } + + faceToSurface[faceI] = hitSurface[0]; + + intersectionFaces.write(faces[faceI], mesh.points()); } } } + +// labelList neiCellSurface(mesh.nFaces()-mesh.nInternalFaces()); +// +// forAll(patches, patchI) +// { +// const polyPatch& pp = patches[patchI]; +// +// if (pp.coupled()) +// { +// forAll(pp, i) +// { +// label faceI = pp.start()+i; +// label ownSurface = cellToSurface[faceOwner[faceI]]; +// neiCellSurface[faceI - mesh.nInternalFaces()] = ownSurface; +// } +// } +// } +// syncTools::swapBoundaryFaceList(mesh, neiCellSurface); +// +// forAll(patches, patchI) +// { +// const polyPatch& pp = patches[patchI]; +// +// if (pp.coupled()) +// { +// forAll(pp, i) +// { +// label faceI = pp.start()+i; +// label ownSurface = cellToSurface[faceOwner[faceI]]; +// label neiSurface = neiCellSurface[faceI-mesh.nInternalFaces()]; +// +// if (faceToSurface[faceI] == -1 && (ownSurface != neiSurface)) +// { +// // Give face the max cell zone +// faceToSurface[faceI] = max(ownSurface, neiSurface); +// } +// } +// } +// } +// // Sync syncTools::syncFaceList(mesh, faceToSurface, maxEqOp