ENH: shm: support for automatic faceZones

This commit is contained in:
Mattijs Janssens
2020-01-16 12:22:02 +00:00
committed by Andrew Heather
parent 0c7140c967
commit b8c2c0acf6
26 changed files with 15695 additions and 322 deletions

View File

@ -597,14 +597,15 @@ private:
) const;
//- Calculate intersections on zoned faces. Return per face -1
// or the index of the surface and the orientation w.r.t. surface
// or the global region of the surface and the orientation
// w.r.t. surface
void getIntersections
(
const labelList& surfacesToTest,
const pointField& neiCc,
const labelList& testFaces,
labelList& namedSurfaceIndex,
labelList& namedSurfaceRegion,
bitSet& posOrientation
) const;
@ -737,13 +738,13 @@ private:
//- Finds zone per cell for cells inside closed named surfaces.
// (uses geometric test for insideness)
// Adapts namedSurfaceIndex so all faces on boundary of cellZone
// Adapts namedSurfaceRegion so all faces on boundary of cellZone
// have corresponding faceZone.
void findCellZoneGeometric
(
const pointField& neiCc,
const labelList& closedNamedSurfaces,
labelList& namedSurfaceIndex,
labelList& namedSurfaceRegion,
const labelList& surfaceToCellZone,
labelList& cellToZone
) const;
@ -781,14 +782,14 @@ private:
//- Finds zone per cell. Uses topological walk with all faces
// marked in unnamedSurfaceRegion (intersections with unnamed
// surfaces) and namedSurfaceIndex (intersections with named
// surfaces) and namedSurfaceRegion (intersections with named
// surfaces) regarded as blocked.
void findCellZoneTopo
(
const label backgroundZoneID,
const pointField& locationsInMesh,
const labelList& unnamedSurfaceRegion,
const labelList& namedSurfaceIndex,
const labelList& namedSurfaceRegion,
const labelList& surfaceToCellZone,
labelList& cellToZone
) const;
@ -800,17 +801,17 @@ private:
const label nErodeCellZones,
const label backgroundZoneID,
const labelList& unnamedSurfaceRegion,
const labelList& namedSurfaceIndex,
const labelList& namedSurfaceRegion,
labelList& cellToZone
) const;
//- Make namedSurfaceIndex consistent with cellToZone
//- Make namedSurfaceRegion consistent with cellToZone
// - clear out any blocked faces inbetween same cell zone.
void makeConsistentFaceIndex
(
const labelList& zoneToNamedSurface,
const labelList& cellToZone,
labelList& namedSurfaceIndex
labelList& namedSurfaceRegion
) const;
//- Calculate cellZone allocation
@ -825,7 +826,7 @@ private:
labelList& cellToZone,
labelList& unnamedRegion1,
labelList& unnamedRegion2,
labelList& namedSurfaceIndex,
labelList& namedSurfaceRegion,
bitSet& posOrientation
) const;

View File

@ -307,7 +307,7 @@ void Foam::meshRefinement::getBafflePatches
labelList cellToZone;
labelList unnamedRegion1;
labelList unnamedRegion2;
labelList namedSurfaceIndex;
labelList namedSurfaceRegion;
{
bitSet posOrientation;
zonify
@ -321,7 +321,7 @@ void Foam::meshRefinement::getBafflePatches
cellToZone,
unnamedRegion1,
unnamedRegion2,
namedSurfaceIndex,
namedSurfaceRegion,
posOrientation
);
}
@ -381,8 +381,8 @@ void Foam::meshRefinement::getBafflePatches
|| (neiZone >= 0 && ownZone != -2)
)
&& (
namedSurfaceIndex.size() == 0
|| namedSurfaceIndex[faceI] == -1
namedSurfaceRegion.size() == 0
|| namedSurfaceRegion[faceI] == -1
)
)
{
@ -426,17 +426,17 @@ Foam::Map<Foam::labelPair> Foam::meshRefinement::getZoneBafflePatches
forAll(surfZones, surfI)
{
const word& faceZoneName = surfZones[surfI].faceZoneName();
const wordList& faceZoneNames = surfZones[surfI].faceZoneNames();
if (faceZoneName.size())
forAll(faceZoneNames, fzi)
{
// Get zone
const word& faceZoneName = faceZoneNames[fzi];
label zoneI = fZones.findZoneID(faceZoneName);
const faceZone& fZone = fZones[zoneI];
// Get patch allocated for zone
label globalRegionI = surfaces_.globalRegion(surfI, 0);
label globalRegionI = surfaces_.globalRegion(surfI, fzi);
labelPair zPatches
(
globalToMasterPatch[globalRegionI],
@ -1370,7 +1370,7 @@ void Foam::meshRefinement::findCellZoneGeometric
(
const pointField& neiCc,
const labelList& closedNamedSurfaces, // indices of closed surfaces
labelList& namedSurfaceIndex, // per face index of named surface
labelList& namedSurfaceRegion, // per face: named surface region
const labelList& surfaceToCellZone, // cell zone index per surface
labelList& cellToZone
@ -1418,11 +1418,9 @@ void Foam::meshRefinement::findCellZoneGeometric
// Count points to test.
label nCandidates = 0;
forAll(namedSurfaceIndex, faceI)
forAll(namedSurfaceRegion, faceI)
{
label surfI = namedSurfaceIndex[faceI];
if (surfI != -1)
if (namedSurfaceRegion[faceI] != -1)
{
if (mesh_.isInternalFace(faceI))
{
@ -1438,11 +1436,9 @@ void Foam::meshRefinement::findCellZoneGeometric
// Collect points.
pointField candidatePoints(nCandidates);
nCandidates = 0;
forAll(namedSurfaceIndex, faceI)
forAll(namedSurfaceRegion, faceI)
{
label surfI = namedSurfaceIndex[faceI];
if (surfI != -1)
if (namedSurfaceRegion[faceI] != -1)
{
label own = faceOwner[faceI];
const point& ownCc = cellCentres[own];
@ -1482,11 +1478,9 @@ void Foam::meshRefinement::findCellZoneGeometric
// 3. Update zone information
nCandidates = 0;
forAll(namedSurfaceIndex, faceI)
forAll(namedSurfaceRegion, faceI)
{
label surfI = namedSurfaceIndex[faceI];
if (surfI != -1)
if (namedSurfaceRegion[faceI] != -1)
{
label own = faceOwner[faceI];
@ -1518,8 +1512,8 @@ void Foam::meshRefinement::findCellZoneGeometric
}
// Adapt the namedSurfaceIndex
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~
// Adapt the namedSurfaceRegion
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// for if any cells were not completely covered.
for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
@ -1527,7 +1521,7 @@ void Foam::meshRefinement::findCellZoneGeometric
label ownZone = cellToZone[mesh_.faceOwner()[faceI]];
label neiZone = cellToZone[mesh_.faceNeighbour()[faceI]];
if (namedSurfaceIndex[faceI] == -1 && (ownZone != neiZone))
if (namedSurfaceRegion[faceI] == -1 && (ownZone != neiZone))
{
// Give face the zone of min cell zone (but only if the
// cellZone originated from a closed, named surface)
@ -1546,12 +1540,14 @@ void Foam::meshRefinement::findCellZoneGeometric
minZone = min(ownZone, neiZone);
}
// Make sure the cellZone originated from a closed surface
// Make sure the cellZone originated from a closed surface. Use
// hardcoded region 0 inside named surface.
label geomSurfI = surfaceToCellZone.find(minZone);
if (geomSurfI != -1)
{
namedSurfaceIndex[faceI] = geomSurfI;
namedSurfaceRegion[faceI] =
surfaces_.globalRegion(geomSurfI, 0);
}
}
}
@ -1573,7 +1569,7 @@ void Foam::meshRefinement::findCellZoneGeometric
label ownZone = cellToZone[mesh_.faceOwner()[faceI]];
label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
if (namedSurfaceIndex[faceI] == -1 && (ownZone != neiZone))
if (namedSurfaceRegion[faceI] == -1 && (ownZone != neiZone))
{
// Give face the min cell zone
label minZone;
@ -1590,12 +1586,14 @@ void Foam::meshRefinement::findCellZoneGeometric
minZone = min(ownZone, neiZone);
}
// Make sure the cellZone originated from a closed surface
// Make sure the cellZone originated from a closed surface.
// Use hardcoded region 0 inside named surface.
label geomSurfI = surfaceToCellZone.find(minZone);
if (geomSurfI != -1)
{
namedSurfaceIndex[faceI] = geomSurfI;
namedSurfaceRegion[faceI] =
surfaces_.globalRegion(geomSurfI, 0);
}
}
}
@ -1603,7 +1601,7 @@ void Foam::meshRefinement::findCellZoneGeometric
}
// Sync
syncTools::syncFaceList(mesh_, namedSurfaceIndex, maxEqOp<label>());
syncTools::syncFaceList(mesh_, namedSurfaceRegion, maxEqOp<label>());
}
@ -1839,7 +1837,7 @@ void Foam::meshRefinement::findCellZoneTopo
const label backgroundZoneID,
const pointField& locationsInMesh,
const labelList& unnamedSurfaceRegion,
const labelList& namedSurfaceIndex,
const labelList& namedSurfaceRegion,
const labelList& surfaceToCellZone,
labelList& cellToZone
) const
@ -1848,7 +1846,7 @@ void Foam::meshRefinement::findCellZoneTopo
// (after all off the unreachable bits of the mesh have been removed).
// This routine splits the mesh into regions, based on the intersection
// with a surface. The problem is that we know the surface which
// (intersected) face belongs to (in namedSurfaceIndex) but we don't
// (intersected) face belongs to (in namedSurfaceRegion) but we don't
// know which side of the face it relates to. So all we are doing here
// is get the correspondence between surface/cellZone and regionSplit
// region. See the logic in calcRegionToZone.
@ -1862,14 +1860,18 @@ void Foam::meshRefinement::findCellZoneTopo
// Assumes:
// - region containing keepPoint does not go into a cellZone
// - all other regions can be found by crossing faces marked in
// namedSurfaceIndex.
// namedSurfaceRegion.
// Analyse regions. Reuse regionsplit
boolList blockedFace(mesh_.nFaces());
forAll(unnamedSurfaceRegion, faceI)
{
if (unnamedSurfaceRegion[faceI] == -1 && namedSurfaceIndex[faceI] == -1)
if
(
unnamedSurfaceRegion[faceI] == -1
&& namedSurfaceRegion[faceI] == -1
)
{
blockedFace[faceI] = false;
}
@ -1878,7 +1880,7 @@ void Foam::meshRefinement::findCellZoneTopo
blockedFace[faceI] = true;
}
}
// No need to sync since namedSurfaceIndex already is synced
// No need to sync since namedSurfaceRegion already is synced
// Set region per cell based on walking
regionSplit cellRegion(mesh_, blockedFace);
@ -1968,13 +1970,16 @@ void Foam::meshRefinement::findCellZoneTopo
for (label faceI = 0; faceI < mesh_.nInternalFaces(); faceI++)
{
label surfI = namedSurfaceIndex[faceI];
label regionI = namedSurfaceRegion[faceI];
// Connected even if no cellZone defined for surface
if (unnamedSurfaceRegion[faceI] == -1 && surfI != -1)
if (unnamedSurfaceRegion[faceI] == -1 && regionI != -1)
{
// Calculate region to zone from cellRegions on either side
// of internal face.
label surfI = surfaces_.whichSurface(regionI).first();
bool changedCell = calcRegionToZone
(
backgroundZoneID,
@ -2008,11 +2013,13 @@ void Foam::meshRefinement::findCellZoneTopo
{
label faceI = pp.start()+i;
label surfI = namedSurfaceIndex[faceI];
label regionI = namedSurfaceRegion[faceI];
// Connected even if no cellZone defined for surface
if (unnamedSurfaceRegion[faceI] == -1 && surfI != -1)
if (unnamedSurfaceRegion[faceI] == -1 && regionI != -1)
{
label surfI = surfaces_.whichSurface(regionI).first();
bool changedCell = calcRegionToZone
(
backgroundZoneID,
@ -2070,7 +2077,7 @@ void Foam::meshRefinement::erodeCellZone
const label nErodeCellZones,
const label backgroundZoneID,
const labelList& unnamedSurfaceRegion,
const labelList& namedSurfaceIndex,
const labelList& namedSurfaceRegion,
labelList& cellToZone
) const
{
@ -2097,7 +2104,7 @@ void Foam::meshRefinement::erodeCellZone
if
(
unnamedSurfaceRegion[facei] == -1
&& namedSurfaceIndex[facei] == -1
&& namedSurfaceRegion[facei] == -1
)
{
label own = mesh_.faceOwner()[facei];
@ -2137,7 +2144,7 @@ void Foam::meshRefinement::erodeCellZone
if
(
unnamedSurfaceRegion[facei] == -1
&& namedSurfaceIndex[facei] == -1
&& namedSurfaceRegion[facei] == -1
)
{
label own = mesh_.faceOwner()[facei];
@ -2174,14 +2181,14 @@ void Foam::meshRefinement::makeConsistentFaceIndex
(
const labelList& surfaceMap,
const labelList& cellToZone,
labelList& namedSurfaceIndex
labelList& namedSurfaceRegion
) const
{
// Make namedSurfaceIndex consistent with cellToZone - clear out any
// Make namedSurfaceRegion consistent with cellToZone - clear out any
// blocked faces inbetween same cell zone (or background (=-1))
// Do not do this for surfaces relating to 'pure' faceZones i.e.
// faceZones without a cellZone. Note that we cannot check here
// for different cellZones on either side but no namedSurfaceIndex
// for different cellZones on either side but no namedSurfaceRegion
// since cellZones can now originate from locationsInMesh as well
// (instead of only through named surfaces)
@ -2192,15 +2199,16 @@ void Foam::meshRefinement::makeConsistentFaceIndex
{
label ownZone = cellToZone[faceOwner[faceI]];
label neiZone = cellToZone[faceNeighbour[faceI]];
label globalI = namedSurfaceRegion[faceI];
if
(
ownZone == neiZone
&& namedSurfaceIndex[faceI] != -1
&& surfaceMap[namedSurfaceIndex[faceI]] == -1
&& globalI != -1
&& surfaceMap[surfaces_.whichSurface(globalI).first()] == -1
)
{
namedSurfaceIndex[faceI] = -1;
namedSurfaceRegion[faceI] = -1;
}
}
@ -2223,15 +2231,16 @@ void Foam::meshRefinement::makeConsistentFaceIndex
label ownZone = cellToZone[faceOwner[faceI]];
label neiZone = neiCellZone[faceI-mesh_.nInternalFaces()];
label globalI = namedSurfaceRegion[faceI];
if
(
ownZone == neiZone
&& namedSurfaceIndex[faceI] != -1
&& surfaceMap[namedSurfaceIndex[faceI]] == -1
&& globalI != -1
&& surfaceMap[surfaces_.whichSurface(globalI).first()] == -1
)
{
namedSurfaceIndex[faceI] = -1;
namedSurfaceRegion[faceI] = -1;
}
}
}
@ -2241,14 +2250,15 @@ void Foam::meshRefinement::makeConsistentFaceIndex
forAll(pp, i)
{
label faceI = pp.start()+i;
label globalI = namedSurfaceRegion[faceI];
if
(
namedSurfaceIndex[faceI] != -1
&& surfaceMap[namedSurfaceIndex[faceI]] == -1
globalI != -1
&& surfaceMap[surfaces_.whichSurface(globalI).first()] == -1
)
{
namedSurfaceIndex[faceI] = -1;
namedSurfaceRegion[faceI] = -1;
}
}
}
@ -2262,12 +2272,12 @@ void Foam::meshRefinement::getIntersections
const pointField& neiCc,
const labelList& testFaces,
labelList& namedSurfaceIndex,
labelList& namedSurfaceRegion,
bitSet& posOrientation
) const
{
namedSurfaceIndex.setSize(mesh_.nFaces());
namedSurfaceIndex = -1;
namedSurfaceRegion.setSize(mesh_.nFaces());
namedSurfaceRegion = -1;
posOrientation.setSize(mesh_.nFaces());
posOrientation = false;
@ -2300,31 +2310,30 @@ void Foam::meshRefinement::getIntersections
// the information already in surfaceIndex_.
labelList surface1;
labelList region1;
List<pointIndexHit> hit1;
vectorField normal1;
labelList surface2;
labelList region2;
List<pointIndexHit> hit2;
vectorField normal2;
{
labelList region1;
labelList region2;
surfaces_.findNearestIntersection
(
surfacesToTest,
start,
end,
surface1,
hit1,
region1,
normal1,
surfaces_.findNearestIntersection
(
surfacesToTest,
start,
end,
surface2,
hit2,
region2,
normal2
);
}
surface1,
hit1,
region1,
normal1,
surface2,
hit2,
region2,
normal2
);
forAll(testFaces, i)
{
@ -2343,20 +2352,32 @@ void Foam::meshRefinement::getIntersections
)
)
{
namedSurfaceIndex[faceI] = surface2[i];
namedSurfaceRegion[faceI] = surfaces_.globalRegion
(
surface2[i],
region2[i]
);
posOrientation.set(faceI, ((area&normal2[i]) > 0));
nSurfFaces[surface2[i]]++;
}
else
{
namedSurfaceIndex[faceI] = surface1[i];
namedSurfaceRegion[faceI] = surfaces_.globalRegion
(
surface1[i],
region1[i]
);
posOrientation.set(faceI, ((area&normal1[i]) > 0));
nSurfFaces[surface1[i]]++;
}
}
else if (surface2[i] != -1)
{
namedSurfaceIndex[faceI] = surface2[i];
namedSurfaceRegion[faceI] = surfaces_.globalRegion
(
surface2[i],
region2[i]
);
posOrientation.set(faceI, ((area&normal2[i]) > 0));
nSurfFaces[surface2[i]]++;
}
@ -2370,7 +2391,7 @@ void Foam::meshRefinement::getIntersections
syncTools::syncFaceList
(
mesh_,
namedSurfaceIndex,
namedSurfaceRegion,
maxEqOp<label>()
);
@ -2399,7 +2420,7 @@ void Foam::meshRefinement::zonify
labelList& cellToZone,
labelList& unnamedRegion1,
labelList& unnamedRegion2,
labelList& namedSurfaceIndex,
labelList& namedSurfaceRegion,
bitSet& posOrientation
) const
{
@ -2408,9 +2429,9 @@ void Foam::meshRefinement::zonify
// -2 : unset
// -1 : not in any zone (zone 'none' or background zone)
// >=0 : zoneID
// namedSurfaceIndex, posOrientation:
// namedSurfaceRegion, posOrientation:
// -1 : face not intersected by named surface
// >=0 : index of named surface
// >=0 : globalRegion (surface+region)
// (and posOrientation: surface normal v.s. face normal)
const PtrList<surfaceZonesInfo>& surfZones = surfaces_.surfZones();
@ -2441,7 +2462,7 @@ void Foam::meshRefinement::zonify
cellToZone.setSize(mesh_.nCells());
cellToZone = -2;
namedSurfaceIndex.clear();
namedSurfaceRegion.clear();
posOrientation.clear();
@ -2466,7 +2487,7 @@ void Foam::meshRefinement::zonify
namedSurfaces,
neiCc,
intersectedFaces(),
namedSurfaceIndex,
namedSurfaceRegion,
posOrientation
);
}
@ -2551,12 +2572,12 @@ void Foam::meshRefinement::zonify
// Stop at unnamed or named surface
labelList allRegion1(mesh_.nFaces(), -1);
forAll(namedSurfaceIndex, faceI)
forAll(namedSurfaceRegion, faceI)
{
allRegion1[faceI] = max
(
unnamedRegion1[faceI],
namedSurfaceIndex[faceI]
namedSurfaceRegion[faceI]
);
}
@ -2595,7 +2616,7 @@ void Foam::meshRefinement::zonify
(
neiCc,
closedNamedSurfaces, // indices of closed surfaces
namedSurfaceIndex, // per face index of named surface
namedSurfaceRegion, // per face index of named surface + region
surfaceToCellZone, // cell zone index per surface
cellToZone
@ -2618,7 +2639,7 @@ void Foam::meshRefinement::zonify
backgroundZoneID,
pointField(0),
unnamedRegion1, // Intersections with unnamed surfaces
namedSurfaceIndex, // Intersections with named surfaces
namedSurfaceRegion, // Intersections with named surfaces
surfaceToCellZone,
cellToZone
);
@ -2635,13 +2656,13 @@ void Foam::meshRefinement::zonify
nErodeCellZones,
backgroundZoneID,
unnamedRegion1,
namedSurfaceIndex,
namedSurfaceRegion,
cellToZone
);
}
// Make sure namedSurfaceIndex is unset inbetween same cell zones.
// Make sure namedSurfaceRegion is unset inbetween same cell zones.
if (!allowFreeStandingZoneFaces)
{
Info<< "Only keeping zone faces inbetween different cellZones."
@ -2668,7 +2689,7 @@ void Foam::meshRefinement::zonify
(
surfaceMap,
cellToZone,
namedSurfaceIndex
namedSurfaceRegion
);
}
}
@ -4476,7 +4497,7 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
// Add any faceZones, cellZones originating from surface to the mesh
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
labelList surfaceToCellZone;
labelList surfaceToFaceZone;
labelListList surfaceToFaceZones;
labelList namedSurfaces(surfaceZonesInfo::getNamedSurfaces(surfZones));
if (namedSurfaces.size())
@ -4487,8 +4508,9 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
label surfI = namedSurfaces[i];
Info<< "Surface : " << surfaces_.names()[surfI] << nl
<< " faceZone : " << surfZones[surfI].faceZoneName() << nl
<< " cellZone : " << surfZones[surfI].cellZoneName() << endl;
<< " faceZones : " << surfZones[surfI].faceZoneNames() << nl
<< " cellZone : " << surfZones[surfI].cellZoneName()
<< endl;
}
Info<< endl;
@ -4499,7 +4521,7 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
namedSurfaces,
mesh_
);
surfaceToFaceZone = surfaceZonesInfo::addFaceZonesToMesh
surfaceToFaceZones = surfaceZonesInfo::addFaceZonesToMesh
(
surfZones,
namedSurfaces,
@ -4516,11 +4538,11 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
// -2 : unset : not allowed!
// -1 : not in any zone (zone 'none')
// >=0: zoneID
// namedSurfaceIndex:
// namedSurfaceRegion:
// -1 : face not intersecting a named surface
// >=0 : index of named surface
labelList cellToZone;
labelList namedSurfaceIndex;
labelList namedSurfaceRegion;
bitSet posOrientation;
{
labelList unnamedRegion1;
@ -4537,24 +4559,26 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::meshRefinement::zonify
cellToZone,
unnamedRegion1,
unnamedRegion2,
namedSurfaceIndex,
namedSurfaceRegion,
posOrientation
);
}
// Convert namedSurfaceIndex (index of named surfaces) to
// Convert namedSurfaceRegion (index of named surfaces) to
// actual faceZone index
//- Per face index of faceZone or -1
labelList faceToZone(mesh_.nFaces(), -1);
forAll(namedSurfaceIndex, faceI)
forAll(namedSurfaceRegion, faceI)
{
label surfI = namedSurfaceIndex[faceI];
if (surfI != -1)
//label surfI = namedSurfaceIndex[faceI];
label globalI = namedSurfaceRegion[faceI];
if (globalI != -1)
{
faceToZone[faceI] = surfaceToFaceZone[surfI];
const labelPair spr = surfaces_.whichSurface(globalI);
faceToZone[faceI] = surfaceToFaceZones[spr.first()][spr.second()];
}
}

View File

@ -123,6 +123,34 @@ Foam::labelList Foam::refinementSurfaces::findHigherLevel
}
Foam::labelList Foam::refinementSurfaces::calcSurfaceIndex
(
const searchableSurfaces& allGeometry,
const labelList& surfaces
)
{
// Determine overall number of global regions
label globalI = 0;
forAll(surfaces, surfI)
{
globalI += allGeometry[surfaces[surfI]].regions().size();
}
labelList regionToSurface(globalI);
globalI = 0;
forAll(surfaces, surfI)
{
const label nLocal = allGeometry[surfaces[surfI]].regions().size();
for (label i = 0; i < nLocal; i++)
{
regionToSurface[globalI++] = surfI;
}
}
return regionToSurface;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::refinementSurfaces::refinementSurfaces
@ -276,7 +304,16 @@ Foam::refinementSurfaces::refinementSurfaces
const searchableSurface& surface = allGeometry_[surfaces_[surfI]];
// Surface zones
surfZones_.set(surfI, new surfaceZonesInfo(surface, dict));
surfZones_.set
(
surfI,
new surfaceZonesInfo
(
surface,
dict,
allGeometry_.regionNames()[surfaces_[surfI]]
)
);
// Global perpendicular angle
if (dict.found("patchInfo"))
@ -433,6 +470,10 @@ Foam::refinementSurfaces::refinementSurfaces
}
// Rework surface specific information into information per global region
regionToSurface_ = calcSurfaceIndex(allGeometry_, surfaces_);
minLevel_.setSize(nRegions);
minLevel_ = 0;
maxLevel_.setSize(nRegions);
@ -539,6 +580,7 @@ Foam::refinementSurfaces::refinementSurfaces
names_(names),
surfZones_(surfZones),
regionOffset_(regionOffset),
regionToSurface_(calcSurfaceIndex(allGeometry, surfaces)),
minLevel_(minLevel),
maxLevel_(maxLevel),
gapLevel_(gapLevel),
@ -558,6 +600,17 @@ Foam::refinementSurfaces::refinementSurfaces
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
Foam::labelPair Foam::refinementSurfaces::whichSurface
(
const label globalRegionI
) const
{
const label surfI = regionToSurface_[globalRegionI];
const label localI = globalRegionI-regionOffset_[surfI];
return labelPair(surfI, localI);
}
// // Count number of triangles per surface region
// Foam::labelList Foam::refinementSurfaces::countRegions(const triSurface& s)
// {

View File

@ -77,9 +77,12 @@ class refinementSurfaces
//- List of surface zone (face and cell zone) information
PtrList<surfaceZonesInfo> surfZones_;
//- From local region number to global region number
//- From surface to starting global region
labelList regionOffset_;
//- From global region number to surface
labelList regionToSurface_;
//- From global region number to refinement level
labelList minLevel_;
@ -121,6 +124,13 @@ class refinementSurfaces
const labelList& surfaceLevel
) const;
//- Calculate global region to surface
static labelList calcSurfaceIndex
(
const searchableSurfaces& allGeometry,
const labelList& surfaces
);
//- No copy construct
refinementSurfaces(const refinementSurfaces&) = delete;
@ -183,7 +193,7 @@ public:
return surfZones_;
}
//- From local region number to global region number
//- From surface to starting global region
const labelList& regionOffset() const
{
return regionOffset_;
@ -252,6 +262,9 @@ public:
return regionOffset_[surfI]+regionI;
}
//- From global region to surface + region
labelPair whichSurface(const label globalRegionI) const;
//- Min level for surface and region on surface
label minLevel(const label surfI, const label regionI) const
{

View File

@ -47,6 +47,18 @@ Foam::surfaceZonesInfo::areaSelectionAlgoNames
});
const Foam::Enum
<
Foam::surfaceZonesInfo::faceZoneNaming
>
Foam::surfaceZonesInfo::faceZoneNamingNames
({
{ faceZoneNaming::NOZONE, "none" },
{ faceZoneNaming::SINGLE, "single" },
{ faceZoneNaming::REGION, "region" }
});
const Foam::Enum
<
Foam::surfaceZonesInfo::faceZoneType
@ -64,18 +76,77 @@ Foam::surfaceZonesInfo::faceZoneTypeNames
Foam::surfaceZonesInfo::surfaceZonesInfo
(
const searchableSurface& surface,
const dictionary& surfacesDict
const dictionary& surfacesDict,
const wordList& regionNames
)
:
faceZoneName_(),
faceZoneNames_(),
cellZoneName_(),
zoneInside_(NONE),
zoneInsidePoint_(point::min),
faceType_(INTERNAL)
{
// Global zone names per surface
if (surfacesDict.readIfPresent("faceZone", faceZoneName_))
const label nRegions = surface.regions().size();
// Old syntax
surfaceZonesInfo::faceZoneNaming namingType = faceZoneNaming::NOZONE;
word namingMethod;
word faceZoneName;
if (surfacesDict.readIfPresent("faceZone", faceZoneName))
{
// Single zone name per surface
if (surfacesDict.found("faceZoneNaming"))
{
FatalIOErrorInFunction(surfacesDict)
<< "Cannot provide both \"faceZone\" and \"faceZoneNaming\""
<< exit(FatalIOError);
}
namingType = faceZoneNaming::SINGLE;
faceZoneNames_.setSize(nRegions, faceZoneName);
}
else if (surfacesDict.readIfPresent("faceZoneNaming", namingMethod))
{
//namingType = faceZoneNamingNames.get("faceZoneNaming", surfacesDict);
namingType = faceZoneNamingNames[namingMethod];
// Generate faceZone names. Maybe make runtime-selection table?
switch (namingType)
{
case faceZoneNaming::NOZONE:
break;
case faceZoneNaming::SINGLE:
{
// Should already be handled above
faceZoneNames_.setSize
(
nRegions,
surfacesDict.get<word>("faceZone")
);
}
break;
case faceZoneNaming::REGION:
{
faceZoneNames_ = regionNames;
}
break;
}
}
if (faceZoneNames_.size())
{
if (faceZoneNames_.size() != nRegions)
{
FatalIOErrorInFunction(surfacesDict)
<< "Number of faceZones (through 'faceZones' keyword)"
<< " does not correspond to the number of regions "
<< nRegions << " in surface " << surface.name()
<< exit(FatalIOError);
}
// Read optional entry to determine inside of faceZone
word method;
@ -116,8 +187,8 @@ Foam::surfaceZonesInfo::surfaceZonesInfo
IOWarningInFunction(surfacesDict)
<< "Illegal entry zoneInside "
<< areaSelectionAlgoNames[zoneInside_]
<< " for faceZone "
<< faceZoneName_
<< " for faceZones "
<< faceZoneNames_
<< " since surface is not closed." << endl;
}
}
@ -125,7 +196,7 @@ Foam::surfaceZonesInfo::surfaceZonesInfo
{
IOWarningInFunction(surfacesDict)
<< "Unused entry zoneInside for faceZone "
<< faceZoneName_
<< faceZoneNames_
<< " since no cellZone specified."
<< endl;
}
@ -142,14 +213,14 @@ Foam::surfaceZonesInfo::surfaceZonesInfo
Foam::surfaceZonesInfo::surfaceZonesInfo
(
const word& faceZoneName,
const wordList& faceZoneNames,
const word& cellZoneName,
const areaSelectionAlgo& zoneInside,
const point& zoneInsidePoint,
const faceZoneType& faceType
)
:
faceZoneName_(faceZoneName),
faceZoneNames_(faceZoneNames),
cellZoneName_(cellZoneName),
zoneInside_(zoneInside),
zoneInsidePoint_(zoneInsidePoint),
@ -159,7 +230,7 @@ Foam::surfaceZonesInfo::surfaceZonesInfo
Foam::surfaceZonesInfo::surfaceZonesInfo(const surfaceZonesInfo& surfZone)
:
faceZoneName_(surfZone.faceZoneName()),
faceZoneNames_(surfZone.faceZoneNames()),
cellZoneName_(surfZone.cellZoneName()),
zoneInside_(surfZone.zoneInside()),
zoneInsidePoint_(surfZone.zoneInsidePoint()),
@ -177,7 +248,7 @@ Foam::labelList Foam::surfaceZonesInfo::getUnnamedSurfaces
label i = 0;
forAll(surfList, surfI)
{
if (surfList[surfI].faceZoneName().empty())
if (surfList[surfI].faceZoneNames().empty())
{
anonymousSurfaces[i++] = surfI;
}
@ -201,7 +272,7 @@ Foam::labelList Foam::surfaceZonesInfo::getNamedSurfaces
if
(
surfList.set(surfI)
&& surfList[surfI].faceZoneName().size()
&& surfList[surfI].faceZoneNames().size()
)
{
namedSurfaces[namedI++] = surfI;
@ -226,7 +297,7 @@ Foam::labelList Foam::surfaceZonesInfo::getStandaloneNamedSurfaces
if
(
surfList.set(surfI)
&& surfList[surfI].faceZoneName().size()
&& surfList[surfI].faceZoneNames().size()
&& !surfList[surfI].cellZoneName().size()
)
{
@ -468,14 +539,14 @@ Foam::label Foam::surfaceZonesInfo::addFaceZone
}
Foam::labelList Foam::surfaceZonesInfo::addFaceZonesToMesh
Foam::labelListList Foam::surfaceZonesInfo::addFaceZonesToMesh
(
const PtrList<surfaceZonesInfo>& surfList,
const labelList& namedSurfaces,
polyMesh& mesh
)
{
labelList surfaceToFaceZone(surfList.size(), -1);
labelListList surfaceToFaceZones(surfList.size());
faceZoneMesh& faceZones = mesh.faceZones();
@ -483,17 +554,23 @@ Foam::labelList Foam::surfaceZonesInfo::addFaceZonesToMesh
{
label surfI = namedSurfaces[i];
const word& faceZoneName = surfList[surfI].faceZoneName();
const wordList& faceZoneNames = surfList[surfI].faceZoneNames();
label zoneI = addFaceZone
(
faceZoneName, //name
labelList(0), //addressing
boolList(0), //flipmap
mesh
);
surfaceToFaceZones[surfI].setSize(faceZoneNames.size(), -1);
forAll(faceZoneNames, j)
{
const word& faceZoneName = faceZoneNames[j];
surfaceToFaceZone[surfI] = zoneI;
label zoneI = addFaceZone
(
faceZoneName, //name
labelList(0), //addressing
boolList(0), //flipmap
mesh
);
surfaceToFaceZones[surfI][j] = zoneI;
}
}
// Check they are synced
@ -515,7 +592,7 @@ Foam::labelList Foam::surfaceZonesInfo::addFaceZonesToMesh
}
}
return surfaceToFaceZone;
return surfaceToFaceZones;
}

View File

@ -73,6 +73,16 @@ public:
static const Enum<areaSelectionAlgo> areaSelectionAlgoNames;
//- How to generate faceZone name
enum faceZoneNaming
{
NOZONE,
SINGLE,
REGION
};
static const Enum<faceZoneNaming> faceZoneNamingNames;
//- What to do with faceZone faces
enum faceZoneType
{
@ -88,8 +98,8 @@ private:
// Private data
//- Per 'interface' surface : name of faceZone to put faces into
word faceZoneName_;
//- Per 'interface' surface : names of faceZones to put faces into
wordList faceZoneNames_;
//- Per 'interface' surface : name of cellZone to put cells into
word cellZoneName_;
@ -117,17 +127,19 @@ public:
// Constructors
//- Construct from surfaces and dictionary
//- Construct from surfaces and dictionary and fully resolved
// region names (for optional automatic faceZone naming)
surfaceZonesInfo
(
const searchableSurface& surface,
const dictionary& surfacesDict
const dictionary& surfacesDict,
const wordList& regionNames
);
//- Construct from components
surfaceZonesInfo
(
const word& faceZoneNames,
const wordList& faceZoneNames,
const word& cellZoneNames,
const areaSelectionAlgo& zoneInside,
const point& zoneInsidePoints,
@ -148,11 +160,11 @@ public:
// Access
//- Per 'interface' surface : empty or name of faceZone to put
// faces into
const word& faceZoneName() const
//- Per 'interface' surface : empty or names of faceZones to put
// faces into (according to region)
const wordList& faceZoneNames() const
{
return faceZoneName_;
return faceZoneNames_;
}
//- Per 'interface' surface : empty or name of cellZone to put
@ -254,7 +266,7 @@ public:
polyMesh& mesh
);
static labelList addFaceZonesToMesh
static labelListList addFaceZonesToMesh
(
const PtrList<surfaceZonesInfo>& surfList,
const labelList& namedSurfaces,

View File

@ -928,11 +928,13 @@ void Foam::snappySnapDriver::preSmoothPatch
// Get (pp-local) indices of points that are both on zone and on patched surface
Foam::labelList Foam::snappySnapDriver::getZoneSurfacePoints
void Foam::snappySnapDriver::getZoneSurfacePoints
(
const fvMesh& mesh,
const indirectPrimitivePatch& pp,
const word& zoneName
const word& zoneName,
bitSet& pointOnZone
)
{
label zonei = mesh.faceZones().findZoneID(zoneName);
@ -950,8 +952,6 @@ Foam::labelList Foam::snappySnapDriver::getZoneSurfacePoints
// Could use PrimitivePatch & localFaces to extract points but might just
// as well do it ourselves.
boolList pointOnZone(pp.nPoints(), false);
forAll(fZone, i)
{
const face& f = mesh.faces()[fZone[i]];
@ -969,8 +969,6 @@ Foam::labelList Foam::snappySnapDriver::getZoneSurfacePoints
}
}
}
return findIndices(pointOnZone, true);
}
@ -1544,124 +1542,128 @@ void Foam::snappySnapDriver::detectNearSurfaces
forAll(zonedSurfaces, i)
{
label zoneSurfi = zonedSurfaces[i];
const word& faceZoneName = surfZones[zoneSurfi].faceZoneName();
const labelList surfacesToTest(1, zoneSurfi);
// Get indices of points both on faceZone and on pp.
labelList zonePointIndices
(
const wordList& faceZoneNames =
surfZones[zoneSurfi].faceZoneNames();
forAll(faceZoneNames, namei)
{
const word& faceZoneName = faceZoneNames[namei];
// Get indices of points both on faceZone and on pp.
bitSet pointOnZone(pp.nPoints());
getZoneSurfacePoints
(
mesh,
pp,
faceZoneName
)
);
faceZoneName,
pointOnZone
);
const labelList zonePointIndices(pointOnZone.toc());
// Do intersection test
labelList surface1;
List<pointIndexHit> hit1;
labelList region1;
vectorField normal1;
// Do intersection test
labelList surface1;
List<pointIndexHit> hit1;
labelList region1;
vectorField normal1;
labelList surface2;
List<pointIndexHit> hit2;
labelList region2;
vectorField normal2;
surfaces.findNearestIntersection
(
surfacesToTest,
pointField(start, zonePointIndices),
pointField(end, zonePointIndices),
labelList surface2;
List<pointIndexHit> hit2;
labelList region2;
vectorField normal2;
surfaces.findNearestIntersection
(
surfacesToTest,
pointField(start, zonePointIndices),
pointField(end, zonePointIndices),
surface1,
hit1,
region1,
normal1,
surface1,
hit1,
region1,
normal1,
surface2,
hit2,
region2,
normal2
);
surface2,
hit2,
region2,
normal2
);
forAll(hit1, i)
{
label pointi = zonePointIndices[i];
// Current location
const point& pt = localPoints[pointi];
bool override = false;
//if (hit1[i].hit())
//{
// if
// (
// meshRefiner_.isGap
// (
// planarCos,
// nearestPoint[pointi],
// nearestNormal[pointi],
// hit1[i].hitPoint(),
// normal1[i]
// )
// )
// {
// disp[pointi] = hit1[i].hitPoint()-pt;
// override = true;
// }
//}
//if (hit2[i].hit())
//{
// if
// (
// meshRefiner_.isGap
// (
// planarCos,
// nearestPoint[pointi],
// nearestNormal[pointi],
// hit2[i].hitPoint(),
// normal2[i]
// )
// )
// {
// disp[pointi] = hit2[i].hitPoint()-pt;
// override = true;
// }
//}
if (hit1[i].hit() && hit2[i].hit())
forAll(hit1, i)
{
if
(
meshRefiner_.isGap
(
planarCos,
hit1[i].hitPoint(),
normal1[i],
hit2[i].hitPoint(),
normal2[i]
)
)
label pointi = zonePointIndices[i];
// Current location
const point& pt = localPoints[pointi];
bool override = false;
//if (hit1[i].hit())
//{
// if
// (
// meshRefiner_.isGap
// (
// planarCos,
// nearestPoint[pointi],
// nearestNormal[pointi],
// hit1[i].hitPoint(),
// normal1[i]
// )
// )
// {
// disp[pointi] = hit1[i].hitPoint()-pt;
// override = true;
// }
//}
//if (hit2[i].hit())
//{
// if
// (
// meshRefiner_.isGap
// (
// planarCos,
// nearestPoint[pointi],
// nearestNormal[pointi],
// hit2[i].hitPoint(),
// normal2[i]
// )
// )
// {
// disp[pointi] = hit2[i].hitPoint()-pt;
// override = true;
// }
//}
if (hit1[i].hit() && hit2[i].hit())
{
if (gapStr.valid())
if
(
meshRefiner_.isGap
(
planarCos,
hit1[i].hitPoint(),
normal1[i],
hit2[i].hitPoint(),
normal2[i]
)
)
{
const point& intPt = hit2[i].hitPoint();
gapStr().write(linePointRef(pt, intPt));
if (gapStr.valid())
{
const point& intPt = hit2[i].hitPoint();
gapStr().write(linePointRef(pt, intPt));
}
disp[pointi] = hit2[i].hitPoint()-pt;
override = true;
}
disp[pointi] = hit2[i].hitPoint()-pt;
override = true;
}
}
if (override && isPatchMasterPoint[pointi])
{
nOverride++;
if (override && isPatchMasterPoint[pointi])
{
nOverride++;
}
}
}
}
@ -1950,15 +1952,15 @@ Foam::vectorField Foam::snappySnapDriver::calcNearestSurface
}
const labelList zonedSurfaces =
surfaceZonesInfo::getNamedSurfaces
(
meshRefiner.surfaces().surfZones()
);
const labelList zonedSurfaces = surfaceZonesInfo::getNamedSurfaces
(
meshRefiner.surfaces().surfZones()
);
// 2. All points on zones to their respective surface
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// (ignoring faceZone subdivision)
// Surfaces with zone information
const PtrList<surfaceZonesInfo>& surfZones = surfaces.surfZones();
@ -1966,26 +1968,27 @@ Foam::vectorField Foam::snappySnapDriver::calcNearestSurface
forAll(zonedSurfaces, i)
{
label surfi = zonedSurfaces[i];
const word& faceZoneName = surfZones[surfi].faceZoneName();
const labelList surfacesToTest(1, surfi);
const label geomi = surfaces.surfaces()[surfi];
const label nRegions =
surfaces.geometry()[geomi].regions().size();
label geomi = surfaces.surfaces()[surfi];
label nRegions = surfaces.geometry()[geomi].regions().size();
const wordList& faceZoneNames =
surfZones[surfi].faceZoneNames();
// Get indices of points both on faceZone and on pp.
labelList zonePointIndices
(
// Get indices of points both on any faceZone and on pp.
bitSet pointOnZone(pp.nPoints());
forAll(faceZoneNames, locali)
{
getZoneSurfacePoints
(
mesh,
pp,
faceZoneName
)
);
faceZoneNames[locali],
pointOnZone
);
}
const labelList zonePointIndices(pointOnZone.toc());
calcNearestSurface
(
@ -2278,9 +2281,12 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::snappySnapDriver::repatchToSurface
forAll(zonedSurfaces, i)
{
const label zoneSurfi = zonedSurfaces[i];
const faceZone& fZone = fZones[surfZones[zoneSurfi].faceZoneName()];
isZonedFace.set(fZone);
const wordList& fZoneNames = surfZones[zoneSurfi].faceZoneNames();
forAll(fZoneNames, i)
{
const faceZone& fZone = fZones[fZoneNames[i]];
isZonedFace.set(fZone);
}
}
}

View File

@ -157,11 +157,13 @@ class snappySnapDriver
labelList getInternalOrBaffleDuplicateFace() const;
//- Get points both on patch and facezone.
static labelList getZoneSurfacePoints
static void getZoneSurfacePoints
(
const fvMesh& mesh,
const indirectPrimitivePatch&,
const word& zoneName
const word& zoneName,
bitSet& pointOnZone
);
//- Get points both on patch and facezone.

View File

@ -263,34 +263,41 @@ void Foam::snappySnapDriver::calcNearestFace
{
label zoneSurfi = zonedSurfaces[i];
const word& faceZoneName = surfZones[zoneSurfi].faceZoneName();
const wordList& faceZoneNames = surfZones[zoneSurfi].faceZoneNames();
// Get indices of faces on pp that are also in zone
label zonei = mesh.faceZones().findZoneID(faceZoneName);
if (zonei == -1)
DynamicList<label> ppFaces;
DynamicList<label> meshFaces;
forAll(faceZoneNames, fzi)
{
FatalErrorInFunction
<< "Problem. Cannot find zone " << faceZoneName
<< exit(FatalError);
}
const faceZone& fZone = mesh.faceZones()[zonei];
const bitSet isZonedFace(mesh.nFaces(), fZone);
DynamicList<label> ppFaces(fZone.size());
DynamicList<label> meshFaces(fZone.size());
forAll(pp.addressing(), i)
{
if (isZonedFace[pp.addressing()[i]])
const word& faceZoneName = faceZoneNames[fzi];
label zonei = mesh.faceZones().findZoneID(faceZoneName);
if (zonei == -1)
{
snapSurf[i] = zoneSurfi;
ppFaces.append(i);
meshFaces.append(pp.addressing()[i]);
FatalErrorInFunction
<< "Problem. Cannot find zone " << faceZoneName
<< exit(FatalError);
}
}
const faceZone& fZone = mesh.faceZones()[zonei];
const bitSet isZonedFace(mesh.nFaces(), fZone);
//Pout<< "For faceZone " << fZone.name()
// << " found " << ppFaces.size() << " out of " << pp.size()
// << endl;
ppFaces.reserve(ppFaces.capacity()+fZone.size());
meshFaces.reserve(meshFaces.capacity()+fZone.size());
forAll(pp.addressing(), i)
{
if (isZonedFace[pp.addressing()[i]])
{
snapSurf[i] = zoneSurfi;
ppFaces.append(i);
meshFaces.append(pp.addressing()[i]);
}
}
//Pout<< "For faceZone " << fZone.name()
// << " found " << ppFaces.size() << " out of " << pp.size()
// << endl;
}
pointField fc
(