zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

ENH: conformalVoronoiMesh: write map from tet mesh to poly

Browse Source
This commit is contained in:
mattijs
2011-07-20 17:43:31 +01:00
parent 68973728a1
commit 4da59cd060
3 changed files with 308 additions and 78 deletions
Download Patch File Download Diff File Expand all files Collapse all files

31
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
View File

@ -684,6 +684,8 @@ private:
labelList& patchStarts,
labelList& procNeighbours,
pointField& cellCentres,
labelList& cellToDelaunayVertex,
labelListList& patchToDelaunayVertex,
bool filterFaces
);
@ -691,6 +693,7 @@ private:
void calcTetMesh
(
pointField& points,
labelList& pointToDelaunayVertex,
faceList& faces,
labelList& owner,
labelList& neighbour,
@ -812,14 +815,15 @@ private:
labelList& patchSizes,
labelList& patchStarts,
labelList& procNeighbours,
labelListList& patchPointPairSlaves,
bool includeEmptyPatches = false
) const;
//- Create the cell centres to use for the mesh
void createCellCentres
(
pointField& cellCentres
) const;
void createCellCentres(pointField& cellCentres) const;
//- Extract all points in vertex-index order
tmp<pointField> allPoints() const;
//- Sort the faces, owner and neighbour lists into
// upper-triangular order. For internal faces only, use
@ -837,6 +841,7 @@ private:
(
List<DynamicList<face> >& patchFaces,
List<DynamicList<label> >& patchOwners,
List<DynamicList<label> >& patchPointPairSlaves,
List<Pair<DynamicList<label> > >& patchSortingIndices
) const;
@ -859,12 +864,12 @@ private:
pointField& pts
) const;
//- Remove dual cells that are not used by any faces
void removeUnusedCells
//- Remove dual cells that are not used by any faces. Return compaction
// map.
labelList removeUnusedCells
(
labelList& owner,
labelList& neighbour,
pointField& cellCentres
labelList& neighbour
) const;
//- Disallow default bitwise copy construct
@ -1008,11 +1013,11 @@ public:
faceList& faces,
labelList& owner,
labelList& neighbour,
wordList& patchTypes,
wordList& patchNames,
labelList& patchSizes,
labelList& patchStarts,
labelList& procNeighbours,
const wordList& patchTypes,
const wordList& patchNames,
const labelList& patchSizes,
const labelList& patchStarts,
const labelList& procNeighbours,
const pointField& cellCentres
) const;

85
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
View File

@ -41,6 +41,8 @@ void Foam::conformalVoronoiMesh::calcDualMesh
labelList& patchStarts,
labelList& procNeighbours,
pointField& cellCentres,
labelList& cellToDelaunayVertex,
labelListList& patchToDelaunayVertex,
bool filterFaces
)
{
@ -274,14 +276,17 @@ void Foam::conformalVoronoiMesh::calcDualMesh
patchSizes,
patchStarts,
procNeighbours,
patchToDelaunayVertex, // from patch face to Delaunay vertex (slavePp)
false
);
// deferredCollapseFaceSet(owner, neighbour, deferredCollapseFaces);
createCellCentres(cellCentres);
cellCentres = allPoints();
removeUnusedCells(owner, neighbour, cellCentres);
cellToDelaunayVertex = removeUnusedCells(owner, neighbour);
cellCentres = pointField(cellCentres, cellToDelaunayVertex);
removeUnusedPoints(faces, points);
@ -292,6 +297,7 @@ void Foam::conformalVoronoiMesh::calcDualMesh
void Foam::conformalVoronoiMesh::calcTetMesh
(
pointField& points,
labelList& pointToDelaunayVertex,
faceList& faces,
labelList& owner,
labelList& neighbour,
@ -306,6 +312,7 @@ void Foam::conformalVoronoiMesh::calcTetMesh
label vertI = 0;
points.setSize(number_of_vertices());
pointToDelaunayVertex.setSize(number_of_vertices());
for
(
@ -318,11 +325,13 @@ void Foam::conformalVoronoiMesh::calcTetMesh
{
vertexMap[vit->index()] = vertI;
points[vertI] = topoint(vit->point());
pointToDelaunayVertex[vertI] = vit->index();
vertI++;
}
}
points.setSize(vertI);
pointToDelaunayVertex.setSize(vertI);
label cellI = 0;
@ -1508,6 +1517,7 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::checkPolyMeshQuality
labelList patchStarts;
labelList procNeighbours;
pointField cellCentres;
labelListList patchToDelaunayVertex;
timeCheck("Start of checkPolyMeshQuality");
@ -1523,12 +1533,15 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::checkPolyMeshQuality
patchSizes,
patchStarts,
procNeighbours,
patchToDelaunayVertex,
false
);
createCellCentres(cellCentres);
//createCellCentres(cellCentres);
cellCentres = allPoints();
removeUnusedCells(owner, neighbour, cellCentres);
labelList cellToDelaunayVertex(removeUnusedCells(owner, neighbour));
cellCentres = pointField(cellCentres, cellToDelaunayVertex);
polyMesh pMesh
(
@ -1887,6 +1900,7 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
labelList& patchSizes,
labelList& patchStarts,
labelList& procNeighbours,
labelListList& patchPointPairSlaves,
bool includeEmptyPatches
) const
{
@ -1966,6 +1980,10 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
List<DynamicList<label> > patchOwners(nPatches, DynamicList<label>(0));
// Per patch face the index of the slave node of the point pair
List<DynamicList<label> > patchPPSlaves(nPatches, DynamicList<label>(0));
faces.setSize(number_of_edges());
owner.setSize(number_of_edges());
@ -2101,6 +2119,16 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
patchFaces[patchIndex].append(newDualFace);
patchOwners[patchIndex].append(own);
// Store the non-internal or boundary point
if (vA->internalOrBoundaryPoint())
{
patchPPSlaves[patchIndex].append(vB->index());
}
else
{
patchPPSlaves[patchIndex].append(vA->index());
}
}
else
{
@ -2139,6 +2167,7 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
(
patchFaces,
patchOwners,
patchPPSlaves,
procPatchSortingIndex
);
@ -2154,6 +2183,13 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
patchFaces,
patchOwners
);
// Return patchPointPairSlaves.setSize(nPatches);
patchPointPairSlaves.setSize(nPatches);
forAll(patchPPSlaves, patchI)
{
patchPointPairSlaves[patchI].transfer(patchPPSlaves[patchI]);
}
}
@ -2162,7 +2198,7 @@ void Foam::conformalVoronoiMesh::createCellCentres
pointField& cellCentres
) const
{
cellCentres.setSize(number_of_vertices());
cellCentres.setSize(number_of_vertices(), point::max);
label vertI = 0;
@ -2185,6 +2221,28 @@ void Foam::conformalVoronoiMesh::createCellCentres
}
Foam::tmp<Foam::pointField> Foam::conformalVoronoiMesh::allPoints() const
{
tmp<pointField> tpts(new pointField(number_of_vertices(), point::max));
pointField& pts = tpts();
for
(
Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
vit != finite_vertices_end();
++vit
)
{
if (vit->internalOrBoundaryPoint())
{
pts[vit->index()] = topoint(vit->point());
}
}
return tpts;
}
void Foam::conformalVoronoiMesh::sortFaces
(
faceList& faces,
@ -2287,6 +2345,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
(
List<DynamicList<face> >& patchFaces,
List<DynamicList<label> >& patchOwners,
List<DynamicList<label> >& patchPointPairSlaves,
List<Pair<DynamicList<label> > >& patchSortingIndices
) const
{
@ -2299,6 +2358,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
{
faceList& faces = patchFaces[patchI];
labelList& owner = patchOwners[patchI];
DynamicList<label>& slaves = patchPointPairSlaves[patchI];
Pair<DynamicList<label> >& sortingIndices = patchSortingIndices[patchI];
@ -2311,6 +2371,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
(
faces.size() != primary.size()
|| owner.size() != primary.size()
|| slaves.size() != primary.size()
)
{
FatalErrorIn
@ -2326,6 +2387,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
<< " for patch " << patchI << nl
<< " faces.size() " << faces.size() << nl
<< " owner.size() " << owner.size() << nl
<< " slaves.size() " << slaves.size() << nl
<< " sortingIndices.first().size() "
<< sortingIndices.first().size()
<< exit(FatalError) << endl;
@ -2344,6 +2406,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
inplaceReorder(oldToNew, secondary);
inplaceReorder(oldToNew, faces);
inplaceReorder(oldToNew, owner);
inplaceReorder(oldToNew, slaves);
// 2) in each block of primary sort by secondary
@ -2401,6 +2464,7 @@ void Foam::conformalVoronoiMesh::sortProcPatches
inplaceReorder(oldToNew, faces);
inplaceReorder(oldToNew, owner);
inplaceReorder(oldToNew, slaves);
}
}
}
@ -2505,11 +2569,10 @@ void Foam::conformalVoronoiMesh::removeUnusedPoints
}
void Foam::conformalVoronoiMesh::removeUnusedCells
Foam::labelList Foam::conformalVoronoiMesh::removeUnusedCells
(
labelList& owner,
labelList& neighbour,
pointField& cellCentres
labelList& neighbour
) const
{
Info<< nl << "Removing unused cells" << endl;
@ -2543,9 +2606,7 @@ void Foam::conformalVoronoiMesh::removeUnusedCells
}
}
inplaceReorder(oldToNew, cellCentres);
cellCentres.setSize(cellI);
labelList newToOld(invert(cellI, oldToNew));
// Find all of the unused cells, create a list of them, then
// subtract one from each owner and neighbour entry for each of
@ -2580,6 +2641,8 @@ void Foam::conformalVoronoiMesh::removeUnusedCells
n -= findLower(unusedCells, n) + 1;
}
}
return newToOld;
}

270
src/mesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C
View File

@ -250,49 +250,12 @@ void Foam::conformalVoronoiMesh::writeMesh
{
writeInternalDelaunayVertices(instance);
if (cvMeshControls().writeTetDualMesh())
{
pointField points;
faceList faces;
labelList owner;
labelList neighbour;
wordList patchTypes;
wordList patchNames;
labelList patchSizes;
labelList patchStarts;
pointField cellCentres;
calcTetMesh
(
points,
faces,
owner,
neighbour,
patchTypes,
patchNames,
patchSizes,
patchStarts
);
labelList procNeighbours(patchNames.size(), -1);
Info<< nl << "Writing tetDualMesh to " << instance << endl;
writeMesh
(
"tetDualMesh",
instance,
points,
faces,
owner,
neighbour,
patchTypes,
patchNames,
patchSizes,
patchStarts,
procNeighbours,
cellCentres
);
}
// Per cell the Delaunay vertex
labelList cellToDelaunayVertex;
// Per patch, per face the Delaunay vertex
labelListList patchToDelaunayVertex;
// Per patch the start of the dual faces
labelList dualPatchStarts;
{
pointField points;
@ -302,7 +265,6 @@ void Foam::conformalVoronoiMesh::writeMesh
wordList patchTypes;
wordList patchNames;
labelList patchSizes;
labelList patchStarts;
labelList procNeighbours;
pointField cellCentres;
@ -315,9 +277,11 @@ void Foam::conformalVoronoiMesh::writeMesh
patchTypes,
patchNames,
patchSizes,
patchStarts,
dualPatchStarts,
procNeighbours,
cellCentres,
cellToDelaunayVertex,
patchToDelaunayVertex,
filterFaces
);
@ -334,6 +298,207 @@ void Foam::conformalVoronoiMesh::writeMesh
patchTypes,
patchNames,
patchSizes,
dualPatchStarts,
procNeighbours,
cellCentres
);
}
if (cvMeshControls().writeTetDualMesh())
{
// Determine map from Delaunay vertex to Dual mesh
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// From all Delaunay vertices to cell (positive index)
// or patch face (negative index)
labelList vertexToDualAddressing(number_of_vertices(), 0);
forAll(cellToDelaunayVertex, cellI)
{
label vertI = cellToDelaunayVertex[cellI];
if (vertexToDualAddressing[vertI] != 0)
{
FatalErrorIn("conformalVoronoiMesh::writeMesh(..)")
<< "Delaunay vertex " << vertI
<< " from cell " << cellI
<< " is already mapped to "
<< vertexToDualAddressing[vertI]
<< exit(FatalError);
}
vertexToDualAddressing[vertI] = cellI+1;
}
forAll(patchToDelaunayVertex, patchI)
{
const labelList& patchVertices = patchToDelaunayVertex[patchI];
forAll(patchVertices, i)
{
label vertI = patchVertices[i];
if (vertexToDualAddressing[vertI] > 0)
{
FatalErrorIn("conformalVoronoiMesh::writeMesh(..)")
<< "Delaunay vertex " << vertI
<< " from patch " << patchI
<< " local index " << i
<< " is already mapped to cell "
<< vertexToDualAddressing[vertI]-1
<< exit(FatalError);
}
// Vertex might be used by multiple faces. Which one to
// use? For now last one wins.
label dualFaceI = dualPatchStarts[patchI]+i;
vertexToDualAddressing[vertI] = -dualFaceI-1;
}
}
// Calculate tet mesh addressing
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
pointField points;
// From tet point back to Delaunay vertex index
labelList pointToDelaunayVertex;
faceList faces;
labelList owner;
labelList neighbour;
wordList patchTypes;
wordList patchNames;
labelList patchSizes;
labelList patchStarts;
pointField cellCentres;
calcTetMesh
(
points,
pointToDelaunayVertex,
faces,
owner,
neighbour,
patchTypes,
patchNames,
patchSizes,
patchStarts
);
// Calculate map from tet points to dual mesh cells/patch faces
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
labelIOList pointDualAddressing
(
IOobject
(
"pointDualAddressing",
instance,
"tetDualMesh"/polyMesh::meshSubDir,
runTime_,
IOobject::NO_READ,
IOobject::AUTO_WRITE,
false
),
UIndirectList<label>
(
vertexToDualAddressing,
pointToDelaunayVertex
)()
);
label pointI = findIndex(pointDualAddressing, -1);
if (pointI != -1)
{
WarningIn
(
"conformalVoronoiMesh::writeMesh\n"
"(\n"
" const fileName& instance,\n"
" bool filterFaces\n"
")\n"
) << "Delaunay vertex " << pointI
<< " does not have a corresponding dual cell." << endl;
}
Info<< "Writing map from tetDualMesh points to Voronoi mesh to "
<< pointDualAddressing.objectPath() << endl;
pointDualAddressing.write();
// Write tet points corresponding to the Voronoi cell/face centre
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
{
// Read Voronoi mesh
fvMesh mesh
(
IOobject
(
Foam::polyMesh::defaultRegion,
instance,
runTime_,
IOobject::MUST_READ
)
);
pointIOField dualPoints
(
IOobject
(
"dualPoints",
instance,
"tetDualMesh"/polyMesh::meshSubDir,
runTime_,
IOobject::NO_READ,
IOobject::AUTO_WRITE,
false
),
points
);
forAll(pointDualAddressing, pointI)
{
label index = pointDualAddressing[pointI];
if (index > 0)
{
label cellI = index-1;
dualPoints[pointI] = mesh.cellCentres()[cellI];
}
else if (index < 0)
{
label faceI = -index-1;
if (faceI >= mesh.nInternalFaces())
{
dualPoints[pointI] = mesh.faceCentres()[faceI];
}
}
}
Info<< "Writing new tetDualMesh points mapped onto Voronoi mesh to "
<< dualPoints.objectPath() << endl
<< "Replace the polyMesh/points with these." << endl;
dualPoints.write();
}
labelList procNeighbours(patchNames.size(), -1);
Info<< nl << "Writing tetDualMesh to " << instance << endl;
writeMesh
(
"tetDualMesh",
instance,
points,
faces,
owner,
neighbour,
patchTypes,
patchNames,
patchSizes,
patchStarts,
procNeighbours,
cellCentres
@ -350,11 +515,11 @@ void Foam::conformalVoronoiMesh::writeMesh
faceList& faces,
labelList& owner,
labelList& neighbour,
wordList& patchTypes,
wordList& patchNames,
labelList& patchSizes,
labelList& patchStarts,
labelList& procNeighbours,
const wordList& patchTypes,
const wordList& patchNames,
const labelList& patchSizes,
const labelList& patchStarts,
const labelList& procNeighbours,
const pointField& cellCentres
) const
{
@ -425,12 +590,9 @@ void Foam::conformalVoronoiMesh::writeMesh
mesh.addFvPatches(patches);
// Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
// mesh.addFvPatches(patches, false);
if (!mesh.write())
{
FatalErrorIn("Foam::conformalVoronoiMesh::writeMesh")
FatalErrorIn("Foam::conformalVoronoiMesh::writeMesh(..)")
<< "Failed writing polyMesh."
<< exit(FatalError);
}