mirror of
https://develop.openfoam.com/Development/openfoam.git
synced 2025-11-28 03:28:01 +00:00
src/sampling/sampledSurface/isoSurface: Correct referencing to temporary fields
Patch provided by Mattijs Janssens Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1487
This commit is contained in:
Henry Weller
@ -37,7 +37,7 @@ License
|
|||||||
|
|
||||||
namespace Foam
|
namespace Foam
|
||||||
{
|
{
|
||||||
defineTypeNameAndDebug(isoSurface, 0);
|
defineTypeNameAndDebug(isoSurface, 0);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@ -58,7 +58,6 @@ bool Foam::isoSurface::noTransform(const tensor& tt) const
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Calculates per face whether couple is collocated.
|
|
||||||
bool Foam::isoSurface::collocatedPatch(const polyPatch& pp)
|
bool Foam::isoSurface::collocatedPatch(const polyPatch& pp)
|
||||||
{
|
{
|
||||||
const coupledPolyPatch& cpp = refCast<const coupledPolyPatch>(pp);
|
const coupledPolyPatch& cpp = refCast<const coupledPolyPatch>(pp);
|
||||||
@ -66,7 +65,6 @@ bool Foam::isoSurface::collocatedPatch(const polyPatch& pp)
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Calculates per face whether couple is collocated.
|
|
||||||
Foam::PackedBoolList Foam::isoSurface::collocatedFaces
|
Foam::PackedBoolList Foam::isoSurface::collocatedFaces
|
||||||
(
|
(
|
||||||
const coupledPolyPatch& pp
|
const coupledPolyPatch& pp
|
||||||
@ -297,7 +295,6 @@ bool Foam::isoSurface::isEdgeOfFaceCut
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Get neighbour value and position.
|
|
||||||
void Foam::isoSurface::getNeighbour
|
void Foam::isoSurface::getNeighbour
|
||||||
(
|
(
|
||||||
const labelList& boundaryRegion,
|
const labelList& boundaryRegion,
|
||||||
@ -331,7 +328,6 @@ void Foam::isoSurface::getNeighbour
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Determine for every face/cell whether it (possibly) generates triangles.
|
|
||||||
void Foam::isoSurface::calcCutTypes
|
void Foam::isoSurface::calcCutTypes
|
||||||
(
|
(
|
||||||
const labelList& boundaryRegion,
|
const labelList& boundaryRegion,
|
||||||
@ -471,43 +467,6 @@ void Foam::isoSurface::calcCutTypes
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Return the two common points between two triangles
|
|
||||||
Foam::labelPair Foam::isoSurface::findCommonPoints
|
|
||||||
(
|
|
||||||
const labelledTri& tri0,
|
|
||||||
const labelledTri& tri1
|
|
||||||
)
|
|
||||||
{
|
|
||||||
labelPair common(-1, -1);
|
|
||||||
|
|
||||||
label fp0 = 0;
|
|
||||||
label fp1 = findIndex(tri1, tri0[fp0]);
|
|
||||||
|
|
||||||
if (fp1 == -1)
|
|
||||||
{
|
|
||||||
fp0 = 1;
|
|
||||||
fp1 = findIndex(tri1, tri0[fp0]);
|
|
||||||
}
|
|
||||||
|
|
||||||
if (fp1 != -1)
|
|
||||||
{
|
|
||||||
// So tri0[fp0] is tri1[fp1]
|
|
||||||
|
|
||||||
// Find next common point
|
|
||||||
label fp0p1 = tri0.fcIndex(fp0);
|
|
||||||
label fp1p1 = tri1.fcIndex(fp1);
|
|
||||||
label fp1m1 = tri1.rcIndex(fp1);
|
|
||||||
|
|
||||||
if (tri0[fp0p1] == tri1[fp1p1] || tri0[fp0p1] == tri1[fp1m1])
|
|
||||||
{
|
|
||||||
common[0] = tri0[fp0];
|
|
||||||
common[1] = tri0[fp0p1];
|
|
||||||
}
|
|
||||||
}
|
|
||||||
return common;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Caculate centre of surface.
|
// Caculate centre of surface.
|
||||||
Foam::point Foam::isoSurface::calcCentre(const triSurface& s)
|
Foam::point Foam::isoSurface::calcCentre(const triSurface& s)
|
||||||
{
|
{
|
||||||
@ -521,75 +480,6 @@ Foam::point Foam::isoSurface::calcCentre(const triSurface& s)
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Replace surface (localPoints, localTris) with single point. Returns
|
|
||||||
// point. Destructs arguments.
|
|
||||||
Foam::pointIndexHit Foam::isoSurface::collapseSurface
|
|
||||||
(
|
|
||||||
pointField& localPoints,
|
|
||||||
DynamicList<labelledTri, 64>& localTris
|
|
||||||
)
|
|
||||||
{
|
|
||||||
pointIndexHit info(false, vector::zero, localTris.size());
|
|
||||||
|
|
||||||
if (localTris.size() == 1)
|
|
||||||
{
|
|
||||||
const labelledTri& tri = localTris[0];
|
|
||||||
info.setPoint(tri.centre(localPoints));
|
|
||||||
info.setHit();
|
|
||||||
}
|
|
||||||
else if (localTris.size() == 2)
|
|
||||||
{
|
|
||||||
// Check if the two triangles share an edge.
|
|
||||||
const labelledTri& tri0 = localTris[0];
|
|
||||||
const labelledTri& tri1 = localTris[0];
|
|
||||||
|
|
||||||
labelPair shared = findCommonPoints(tri0, tri1);
|
|
||||||
|
|
||||||
if (shared[0] != -1)
|
|
||||||
{
|
|
||||||
info.setPoint
|
|
||||||
(
|
|
||||||
0.5
|
|
||||||
* (
|
|
||||||
tri0.centre(localPoints)
|
|
||||||
+ tri1.centre(localPoints)
|
|
||||||
)
|
|
||||||
);
|
|
||||||
info.setHit();
|
|
||||||
}
|
|
||||||
}
|
|
||||||
else if (localTris.size())
|
|
||||||
{
|
|
||||||
// Check if single region. Rare situation.
|
|
||||||
triSurface surf
|
|
||||||
(
|
|
||||||
localTris,
|
|
||||||
geometricSurfacePatchList(0),
|
|
||||||
localPoints,
|
|
||||||
true
|
|
||||||
);
|
|
||||||
localTris.clearStorage();
|
|
||||||
|
|
||||||
labelList faceZone;
|
|
||||||
label nZones = surf.markZones
|
|
||||||
(
|
|
||||||
boolList(surf.nEdges(), false),
|
|
||||||
faceZone
|
|
||||||
);
|
|
||||||
|
|
||||||
if (nZones == 1)
|
|
||||||
{
|
|
||||||
info.setPoint(calcCentre(surf));
|
|
||||||
info.setHit();
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
return info;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Determine per cell centre whether all the intersections get collapsed
|
|
||||||
// to a single point
|
|
||||||
void Foam::isoSurface::calcSnappedCc
|
void Foam::isoSurface::calcSnappedCc
|
||||||
(
|
(
|
||||||
const labelList& boundaryRegion,
|
const labelList& boundaryRegion,
|
||||||
@ -755,8 +645,6 @@ void Foam::isoSurface::calcSnappedCc
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Determine per meshpoint whether all the intersections get collapsed
|
|
||||||
// to a single point
|
|
||||||
void Foam::isoSurface::calcSnappedPoint
|
void Foam::isoSurface::calcSnappedPoint
|
||||||
(
|
(
|
||||||
const PackedBoolList& isBoundaryPoint,
|
const PackedBoolList& isBoundaryPoint,
|
||||||
@ -1131,7 +1019,6 @@ Foam::triSurface Foam::isoSurface::stitchTriPoints
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Does face use valid vertices?
|
|
||||||
bool Foam::isoSurface::validTri(const triSurface& surf, const label faceI)
|
bool Foam::isoSurface::validTri(const triSurface& surf, const label faceI)
|
||||||
{
|
{
|
||||||
// Simple check on indices ok.
|
// Simple check on indices ok.
|
||||||
@ -1202,435 +1089,6 @@ bool Foam::isoSurface::validTri(const triSurface& surf, const label faceI)
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
void Foam::isoSurface::calcAddressing
|
|
||||||
(
|
|
||||||
const triSurface& surf,
|
|
||||||
List<FixedList<label, 3>>& faceEdges,
|
|
||||||
labelList& edgeFace0,
|
|
||||||
labelList& edgeFace1,
|
|
||||||
Map<labelList>& edgeFacesRest
|
|
||||||
) const
|
|
||||||
{
|
|
||||||
const pointField& points = surf.points();
|
|
||||||
|
|
||||||
pointField edgeCentres(3*surf.size());
|
|
||||||
label edgeI = 0;
|
|
||||||
forAll(surf, triI)
|
|
||||||
{
|
|
||||||
const labelledTri& tri = surf[triI];
|
|
||||||
edgeCentres[edgeI++] = 0.5*(points[tri[0]]+points[tri[1]]);
|
|
||||||
edgeCentres[edgeI++] = 0.5*(points[tri[1]]+points[tri[2]]);
|
|
||||||
edgeCentres[edgeI++] = 0.5*(points[tri[2]]+points[tri[0]]);
|
|
||||||
}
|
|
||||||
|
|
||||||
pointField mergedCentres;
|
|
||||||
labelList oldToMerged;
|
|
||||||
bool hasMerged = mergePoints
|
|
||||||
(
|
|
||||||
edgeCentres,
|
|
||||||
mergeDistance_,
|
|
||||||
false,
|
|
||||||
oldToMerged,
|
|
||||||
mergedCentres
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
|
||||||
{
|
|
||||||
Pout<< "isoSurface : detected "
|
|
||||||
<< mergedCentres.size()
|
|
||||||
<< " geometric edges on " << surf.size() << " triangles." << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (!hasMerged)
|
|
||||||
{
|
|
||||||
return;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Determine faceEdges
|
|
||||||
faceEdges.setSize(surf.size());
|
|
||||||
edgeI = 0;
|
|
||||||
forAll(surf, triI)
|
|
||||||
{
|
|
||||||
faceEdges[triI][0] = oldToMerged[edgeI++];
|
|
||||||
faceEdges[triI][1] = oldToMerged[edgeI++];
|
|
||||||
faceEdges[triI][2] = oldToMerged[edgeI++];
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Determine edgeFaces
|
|
||||||
edgeFace0.setSize(mergedCentres.size());
|
|
||||||
edgeFace0 = -1;
|
|
||||||
edgeFace1.setSize(mergedCentres.size());
|
|
||||||
edgeFace1 = -1;
|
|
||||||
edgeFacesRest.clear();
|
|
||||||
|
|
||||||
// Overflow edge faces for geometric shared edges that turned
|
|
||||||
// out to be different anyway.
|
|
||||||
EdgeMap<labelList> extraEdgeFaces(mergedCentres.size()/100);
|
|
||||||
|
|
||||||
forAll(oldToMerged, oldEdgeI)
|
|
||||||
{
|
|
||||||
label triI = oldEdgeI / 3;
|
|
||||||
label edgeI = oldToMerged[oldEdgeI];
|
|
||||||
|
|
||||||
if (edgeFace0[edgeI] == -1)
|
|
||||||
{
|
|
||||||
// First triangle for edge
|
|
||||||
edgeFace0[edgeI] = triI;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
//- Check that the two triangles actually topologically
|
|
||||||
// share an edge
|
|
||||||
const labelledTri& prevTri = surf[edgeFace0[edgeI]];
|
|
||||||
const labelledTri& tri = surf[triI];
|
|
||||||
|
|
||||||
label fp = oldEdgeI % 3;
|
|
||||||
|
|
||||||
edge e(tri[fp], tri[tri.fcIndex(fp)]);
|
|
||||||
|
|
||||||
label prevTriIndex = -1;
|
|
||||||
|
|
||||||
forAll(prevTri, i)
|
|
||||||
{
|
|
||||||
if (edge(prevTri[i], prevTri[prevTri.fcIndex(i)]) == e)
|
|
||||||
{
|
|
||||||
prevTriIndex = i;
|
|
||||||
break;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
if (prevTriIndex == -1)
|
|
||||||
{
|
|
||||||
// Different edge. Store for later.
|
|
||||||
EdgeMap<labelList>::iterator iter = extraEdgeFaces.find(e);
|
|
||||||
|
|
||||||
if (iter != extraEdgeFaces.end())
|
|
||||||
{
|
|
||||||
labelList& eFaces = iter();
|
|
||||||
label sz = eFaces.size();
|
|
||||||
eFaces.setSize(sz+1);
|
|
||||||
eFaces[sz] = triI;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
extraEdgeFaces.insert(e, labelList(1, triI));
|
|
||||||
}
|
|
||||||
}
|
|
||||||
else if (edgeFace1[edgeI] == -1)
|
|
||||||
{
|
|
||||||
edgeFace1[edgeI] = triI;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
//WarningInFunction
|
|
||||||
// << "Edge " << edgeI << " with centre "
|
|
||||||
// << mergedCentres[edgeI]
|
|
||||||
// << " used by more than two triangles: "
|
|
||||||
// << edgeFace0[edgeI] << ", "
|
|
||||||
// << edgeFace1[edgeI] << " and " << triI << endl;
|
|
||||||
Map<labelList>::iterator iter = edgeFacesRest.find(edgeI);
|
|
||||||
|
|
||||||
if (iter != edgeFacesRest.end())
|
|
||||||
{
|
|
||||||
labelList& eFaces = iter();
|
|
||||||
label sz = eFaces.size();
|
|
||||||
eFaces.setSize(sz+1);
|
|
||||||
eFaces[sz] = triI;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
edgeFacesRest.insert(edgeI, labelList(1, triI));
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
// Add extraEdgeFaces
|
|
||||||
edgeI = edgeFace0.size();
|
|
||||||
|
|
||||||
edgeFace0.setSize(edgeI + extraEdgeFaces.size());
|
|
||||||
edgeFace1.setSize(edgeI + extraEdgeFaces.size(), -1);
|
|
||||||
|
|
||||||
forAllConstIter(EdgeMap<labelList>, extraEdgeFaces, iter)
|
|
||||||
{
|
|
||||||
const labelList& eFaces = iter();
|
|
||||||
|
|
||||||
// The current edge will become edgeI. Replace all occurrences in
|
|
||||||
// faceEdges
|
|
||||||
forAll(eFaces, i)
|
|
||||||
{
|
|
||||||
label triI = eFaces[i];
|
|
||||||
const labelledTri& tri = surf[triI];
|
|
||||||
|
|
||||||
FixedList<label, 3>& fEdges = faceEdges[triI];
|
|
||||||
forAll(tri, fp)
|
|
||||||
{
|
|
||||||
edge e(tri[fp], tri[tri.fcIndex(fp)]);
|
|
||||||
|
|
||||||
if (e == iter.key())
|
|
||||||
{
|
|
||||||
fEdges[fp] = edgeI;
|
|
||||||
break;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Add face to edgeFaces
|
|
||||||
|
|
||||||
edgeFace0[edgeI] = eFaces[0];
|
|
||||||
|
|
||||||
if (eFaces.size() >= 2)
|
|
||||||
{
|
|
||||||
edgeFace1[edgeI] = eFaces[1];
|
|
||||||
|
|
||||||
if (eFaces.size() > 2)
|
|
||||||
{
|
|
||||||
edgeFacesRest.insert
|
|
||||||
(
|
|
||||||
edgeI,
|
|
||||||
SubList<label>(eFaces, eFaces.size()-2, 2)
|
|
||||||
);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
edgeI++;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
void Foam::isoSurface::walkOrientation
|
|
||||||
(
|
|
||||||
const triSurface& surf,
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const label seedTriI,
|
|
||||||
labelList& flipState
|
|
||||||
)
|
|
||||||
{
|
|
||||||
// Do walk for consistent orientation.
|
|
||||||
DynamicList<label> changedFaces(surf.size());
|
|
||||||
|
|
||||||
changedFaces.append(seedTriI);
|
|
||||||
|
|
||||||
while (changedFaces.size())
|
|
||||||
{
|
|
||||||
DynamicList<label> newChangedFaces(changedFaces.size());
|
|
||||||
|
|
||||||
forAll(changedFaces, i)
|
|
||||||
{
|
|
||||||
label triI = changedFaces[i];
|
|
||||||
const labelledTri& tri = surf[triI];
|
|
||||||
const FixedList<label, 3>& fEdges = faceEdges[triI];
|
|
||||||
|
|
||||||
forAll(fEdges, fp)
|
|
||||||
{
|
|
||||||
label edgeI = fEdges[fp];
|
|
||||||
|
|
||||||
// my points:
|
|
||||||
label p0 = tri[fp];
|
|
||||||
label p1 = tri[tri.fcIndex(fp)];
|
|
||||||
|
|
||||||
label nbrI =
|
|
||||||
(
|
|
||||||
edgeFace0[edgeI] != triI
|
|
||||||
? edgeFace0[edgeI]
|
|
||||||
: edgeFace1[edgeI]
|
|
||||||
);
|
|
||||||
|
|
||||||
if (nbrI != -1 && flipState[nbrI] == -1)
|
|
||||||
{
|
|
||||||
const labelledTri& nbrTri = surf[nbrI];
|
|
||||||
|
|
||||||
// nbr points
|
|
||||||
label nbrFp = findIndex(nbrTri, p0);
|
|
||||||
|
|
||||||
if (nbrFp == -1)
|
|
||||||
{
|
|
||||||
FatalErrorInFunction
|
|
||||||
<< "triI:" << triI
|
|
||||||
<< " tri:" << tri
|
|
||||||
<< " p0:" << p0
|
|
||||||
<< " p1:" << p1
|
|
||||||
<< " fEdges:" << fEdges
|
|
||||||
<< " edgeI:" << edgeI
|
|
||||||
<< " edgeFace0:" << edgeFace0[edgeI]
|
|
||||||
<< " edgeFace1:" << edgeFace1[edgeI]
|
|
||||||
<< " nbrI:" << nbrI
|
|
||||||
<< " nbrTri:" << nbrTri
|
|
||||||
<< abort(FatalError);
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
label nbrP1 = nbrTri[nbrTri.rcIndex(nbrFp)];
|
|
||||||
|
|
||||||
bool sameOrientation = (p1 == nbrP1);
|
|
||||||
|
|
||||||
if (flipState[triI] == 0)
|
|
||||||
{
|
|
||||||
flipState[nbrI] = (sameOrientation ? 0 : 1);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
flipState[nbrI] = (sameOrientation ? 1 : 0);
|
|
||||||
}
|
|
||||||
newChangedFaces.append(nbrI);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
changedFaces.transfer(newChangedFaces);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
void Foam::isoSurface::orientSurface
|
|
||||||
(
|
|
||||||
triSurface& surf,
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const Map<labelList>& edgeFacesRest
|
|
||||||
)
|
|
||||||
{
|
|
||||||
// -1 : unvisited
|
|
||||||
// 0 : leave as is
|
|
||||||
// 1 : flip
|
|
||||||
labelList flipState(surf.size(), -1);
|
|
||||||
|
|
||||||
label seedTriI = 0;
|
|
||||||
|
|
||||||
while (true)
|
|
||||||
{
|
|
||||||
// Find first unvisited triangle
|
|
||||||
for
|
|
||||||
(
|
|
||||||
;
|
|
||||||
seedTriI < surf.size() && flipState[seedTriI] != -1;
|
|
||||||
seedTriI++
|
|
||||||
)
|
|
||||||
{}
|
|
||||||
|
|
||||||
if (seedTriI == surf.size())
|
|
||||||
{
|
|
||||||
break;
|
|
||||||
}
|
|
||||||
|
|
||||||
// Note: Determine orientation of seedTriI?
|
|
||||||
// for now assume it is ok
|
|
||||||
flipState[seedTriI] = 0;
|
|
||||||
|
|
||||||
walkOrientation
|
|
||||||
(
|
|
||||||
surf,
|
|
||||||
faceEdges,
|
|
||||||
edgeFace0,
|
|
||||||
edgeFace1,
|
|
||||||
seedTriI,
|
|
||||||
flipState
|
|
||||||
);
|
|
||||||
}
|
|
||||||
|
|
||||||
// Do actual flipping
|
|
||||||
surf.clearOut();
|
|
||||||
forAll(surf, triI)
|
|
||||||
{
|
|
||||||
if (flipState[triI] == 1)
|
|
||||||
{
|
|
||||||
labelledTri tri(surf[triI]);
|
|
||||||
|
|
||||||
surf[triI][0] = tri[0];
|
|
||||||
surf[triI][1] = tri[2];
|
|
||||||
surf[triI][2] = tri[1];
|
|
||||||
}
|
|
||||||
else if (flipState[triI] == -1)
|
|
||||||
{
|
|
||||||
FatalErrorInFunction
|
|
||||||
<< "problem" << abort(FatalError);
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Checks if triangle is connected through edgeI only.
|
|
||||||
bool Foam::isoSurface::danglingTriangle
|
|
||||||
(
|
|
||||||
const FixedList<label, 3>& fEdges,
|
|
||||||
const labelList& edgeFace1
|
|
||||||
)
|
|
||||||
{
|
|
||||||
label nOpen = 0;
|
|
||||||
forAll(fEdges, i)
|
|
||||||
{
|
|
||||||
if (edgeFace1[fEdges[i]] == -1)
|
|
||||||
{
|
|
||||||
nOpen++;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
if (nOpen == 1 || nOpen == 2 || nOpen == 3)
|
|
||||||
{
|
|
||||||
return true;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
return false;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
// Mark triangles to keep. Returns number of dangling triangles.
|
|
||||||
Foam::label Foam::isoSurface::markDanglingTriangles
|
|
||||||
(
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const Map<labelList>& edgeFacesRest,
|
|
||||||
boolList& keepTriangles
|
|
||||||
)
|
|
||||||
{
|
|
||||||
keepTriangles.setSize(faceEdges.size());
|
|
||||||
keepTriangles = true;
|
|
||||||
|
|
||||||
label nDangling = 0;
|
|
||||||
|
|
||||||
// Remove any dangling triangles
|
|
||||||
forAllConstIter(Map<labelList>, edgeFacesRest, iter)
|
|
||||||
{
|
|
||||||
// These are all the non-manifold edges. Filter out all triangles
|
|
||||||
// with only one connected edge (= this edge)
|
|
||||||
|
|
||||||
label edgeI = iter.key();
|
|
||||||
const labelList& otherEdgeFaces = iter();
|
|
||||||
|
|
||||||
// Remove all dangling triangles
|
|
||||||
if (danglingTriangle(faceEdges[edgeFace0[edgeI]], edgeFace1))
|
|
||||||
{
|
|
||||||
keepTriangles[edgeFace0[edgeI]] = false;
|
|
||||||
nDangling++;
|
|
||||||
}
|
|
||||||
if (danglingTriangle(faceEdges[edgeFace1[edgeI]], edgeFace1))
|
|
||||||
{
|
|
||||||
keepTriangles[edgeFace1[edgeI]] = false;
|
|
||||||
nDangling++;
|
|
||||||
}
|
|
||||||
forAll(otherEdgeFaces, i)
|
|
||||||
{
|
|
||||||
label triI = otherEdgeFaces[i];
|
|
||||||
if (danglingTriangle(faceEdges[triI], edgeFace1))
|
|
||||||
{
|
|
||||||
keepTriangles[triI] = false;
|
|
||||||
nDangling++;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
return nDangling;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
Foam::triSurface Foam::isoSurface::subsetMesh
|
Foam::triSurface Foam::isoSurface::subsetMesh
|
||||||
(
|
(
|
||||||
const triSurface& s,
|
const triSurface& s,
|
||||||
@ -1953,57 +1411,58 @@ Foam::isoSurface::isoSurface
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
DynamicList<point> triPoints(nCutCells_);
|
|
||||||
DynamicList<label> triMeshCells(nCutCells_);
|
|
||||||
|
|
||||||
generateTriPoints
|
|
||||||
(
|
|
||||||
cValsPtr_(),
|
|
||||||
pVals_,
|
|
||||||
|
|
||||||
meshC,
|
|
||||||
mesh_.points(),
|
|
||||||
|
|
||||||
snappedPoints,
|
|
||||||
snappedCc,
|
|
||||||
snappedPoint,
|
|
||||||
|
|
||||||
triPoints,
|
|
||||||
triMeshCells
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
|
||||||
{
|
{
|
||||||
Pout<< "isoSurface : generated " << triMeshCells.size()
|
DynamicList<point> triPoints(3*nCutCells_);
|
||||||
<< " unmerged triangles from " << triPoints.size()
|
DynamicList<label> triMeshCells(nCutCells_);
|
||||||
<< " unmerged points." << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
generateTriPoints
|
||||||
// Merge points and compact out non-valid triangles
|
|
||||||
labelList triMap; // merged to unmerged triangle
|
|
||||||
triSurface::operator=
|
|
||||||
(
|
|
||||||
stitchTriPoints
|
|
||||||
(
|
(
|
||||||
true, // check for duplicate tris
|
cValsPtr_(),
|
||||||
triPoints,
|
pVals_,
|
||||||
triPointMergeMap_, // unmerged to merged point
|
|
||||||
triMap
|
|
||||||
)
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
meshC,
|
||||||
{
|
mesh_.points(),
|
||||||
Pout<< "isoSurface : generated " << triMap.size()
|
|
||||||
<< " merged triangles." << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
meshCells_.setSize(triMap.size());
|
snappedPoints,
|
||||||
forAll(triMap, i)
|
snappedCc,
|
||||||
{
|
snappedPoint,
|
||||||
meshCells_[i] = triMeshCells[triMap[i]];
|
|
||||||
|
triPoints, // 3 points of the triangle
|
||||||
|
triMeshCells // per triangle the originating cell
|
||||||
|
);
|
||||||
|
|
||||||
|
if (debug)
|
||||||
|
{
|
||||||
|
Pout<< "isoSurface : generated " << triMeshCells.size()
|
||||||
|
<< " unmerged triangles from " << triPoints.size()
|
||||||
|
<< " unmerged points." << endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
// Merge points and compact out non-valid triangles
|
||||||
|
labelList triMap; // merged to unmerged triangle
|
||||||
|
triSurface::operator=
|
||||||
|
(
|
||||||
|
stitchTriPoints
|
||||||
|
(
|
||||||
|
true, // check for duplicate tris
|
||||||
|
triPoints,
|
||||||
|
triPointMergeMap_, // unmerged to merged point
|
||||||
|
triMap
|
||||||
|
)
|
||||||
|
);
|
||||||
|
|
||||||
|
if (debug)
|
||||||
|
{
|
||||||
|
Pout<< "isoSurface : generated " << triMap.size()
|
||||||
|
<< " merged triangles." << endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
meshCells_.setSize(triMap.size());
|
||||||
|
forAll(triMap, i)
|
||||||
|
{
|
||||||
|
meshCells_[i] = triMeshCells[triMap[i]];
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
if (debug)
|
if (debug)
|
||||||
@ -2019,70 +1478,6 @@ Foam::isoSurface::isoSurface
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
if (false)
|
|
||||||
{
|
|
||||||
List<FixedList<label, 3>> faceEdges;
|
|
||||||
labelList edgeFace0, edgeFace1;
|
|
||||||
Map<labelList> edgeFacesRest;
|
|
||||||
|
|
||||||
|
|
||||||
while (true)
|
|
||||||
{
|
|
||||||
// Calculate addressing
|
|
||||||
calcAddressing
|
|
||||||
(
|
|
||||||
*this,
|
|
||||||
faceEdges,
|
|
||||||
edgeFace0,
|
|
||||||
edgeFace1,
|
|
||||||
edgeFacesRest
|
|
||||||
);
|
|
||||||
|
|
||||||
// See if any dangling triangles
|
|
||||||
boolList keepTriangles;
|
|
||||||
label nDangling = markDanglingTriangles
|
|
||||||
(
|
|
||||||
faceEdges,
|
|
||||||
edgeFace0,
|
|
||||||
edgeFace1,
|
|
||||||
edgeFacesRest,
|
|
||||||
keepTriangles
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
|
||||||
{
|
|
||||||
Pout<< "isoSurface : detected " << nDangling
|
|
||||||
<< " dangling triangles." << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
if (nDangling == 0)
|
|
||||||
{
|
|
||||||
break;
|
|
||||||
}
|
|
||||||
|
|
||||||
// Create face map (new to old)
|
|
||||||
labelList subsetTriMap(findIndices(keepTriangles, true));
|
|
||||||
|
|
||||||
labelList subsetPointMap;
|
|
||||||
labelList reversePointMap;
|
|
||||||
triSurface::operator=
|
|
||||||
(
|
|
||||||
subsetMesh
|
|
||||||
(
|
|
||||||
*this,
|
|
||||||
subsetTriMap,
|
|
||||||
reversePointMap,
|
|
||||||
subsetPointMap
|
|
||||||
)
|
|
||||||
);
|
|
||||||
meshCells_ = labelField(meshCells_, subsetTriMap);
|
|
||||||
inplaceRenumber(reversePointMap, triPointMergeMap_);
|
|
||||||
}
|
|
||||||
|
|
||||||
orientSurface(*this, faceEdges, edgeFace0, edgeFace1, edgeFacesRest);
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
if (debug)
|
if (debug)
|
||||||
{
|
{
|
||||||
fileName stlFile = mesh_.time().path() + ".stl";
|
fileName stlFile = mesh_.time().path() + ".stl";
|
||||||
|
|||||||
@ -173,6 +173,7 @@ class isoSurface
|
|||||||
const bool neiLower
|
const bool neiLower
|
||||||
) const;
|
) const;
|
||||||
|
|
||||||
|
//- Get neighbour value and position.
|
||||||
void getNeighbour
|
void getNeighbour
|
||||||
(
|
(
|
||||||
const labelList& boundaryRegion,
|
const labelList& boundaryRegion,
|
||||||
@ -184,7 +185,8 @@ class isoSurface
|
|||||||
point& nbrPoint
|
point& nbrPoint
|
||||||
) const;
|
) const;
|
||||||
|
|
||||||
//- Set faceCutType,cellCutType.
|
//- Determine for every face/cell whether it (possibly) generates
|
||||||
|
// triangles.
|
||||||
void calcCutTypes
|
void calcCutTypes
|
||||||
(
|
(
|
||||||
const labelList& boundaryRegion,
|
const labelList& boundaryRegion,
|
||||||
@ -193,20 +195,8 @@ class isoSurface
|
|||||||
const scalarField& pVals
|
const scalarField& pVals
|
||||||
);
|
);
|
||||||
|
|
||||||
static labelPair findCommonPoints
|
|
||||||
(
|
|
||||||
const labelledTri&,
|
|
||||||
const labelledTri&
|
|
||||||
);
|
|
||||||
|
|
||||||
static point calcCentre(const triSurface&);
|
static point calcCentre(const triSurface&);
|
||||||
|
|
||||||
static pointIndexHit collapseSurface
|
|
||||||
(
|
|
||||||
pointField& localPoints,
|
|
||||||
DynamicList<labelledTri, 64>& localTris
|
|
||||||
);
|
|
||||||
|
|
||||||
//- Determine per cc whether all near cuts can be snapped to single
|
//- Determine per cc whether all near cuts can be snapped to single
|
||||||
// point.
|
// point.
|
||||||
void calcSnappedCc
|
void calcSnappedCc
|
||||||
@ -257,6 +247,9 @@ class isoSurface
|
|||||||
const Type& snapP1
|
const Type& snapP1
|
||||||
) const;
|
) const;
|
||||||
|
|
||||||
|
|
||||||
|
//- Note: cannot use simpler isoSurfaceCell::generateTriPoints since
|
||||||
|
// the need here to sometimes pass in remote 'snappedPoints'
|
||||||
template<class Type>
|
template<class Type>
|
||||||
void generateTriPoints
|
void generateTriPoints
|
||||||
(
|
(
|
||||||
@ -323,6 +316,14 @@ class isoSurface
|
|||||||
DynamicList<label>& triMeshCells
|
DynamicList<label>& triMeshCells
|
||||||
) const;
|
) const;
|
||||||
|
|
||||||
|
template<class Type>
|
||||||
|
static tmp<Field<Type>> interpolate
|
||||||
|
(
|
||||||
|
const label nPoints,
|
||||||
|
const labelList& triPointMergeMap,
|
||||||
|
const DynamicList<Type>& unmergedValues
|
||||||
|
);
|
||||||
|
|
||||||
triSurface stitchTriPoints
|
triSurface stitchTriPoints
|
||||||
(
|
(
|
||||||
const bool checkDuplicates,
|
const bool checkDuplicates,
|
||||||
@ -334,54 +335,6 @@ class isoSurface
|
|||||||
//- Check single triangle for (topological) validity
|
//- Check single triangle for (topological) validity
|
||||||
static bool validTri(const triSurface&, const label);
|
static bool validTri(const triSurface&, const label);
|
||||||
|
|
||||||
//- Determine edge-face addressing
|
|
||||||
void calcAddressing
|
|
||||||
(
|
|
||||||
const triSurface& surf,
|
|
||||||
List<FixedList<label, 3>>& faceEdges,
|
|
||||||
labelList& edgeFace0,
|
|
||||||
labelList& edgeFace1,
|
|
||||||
Map<labelList>& edgeFacesRest
|
|
||||||
) const;
|
|
||||||
|
|
||||||
//- Determine orientation
|
|
||||||
static void walkOrientation
|
|
||||||
(
|
|
||||||
const triSurface& surf,
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const label seedTriI,
|
|
||||||
labelList& flipState
|
|
||||||
);
|
|
||||||
|
|
||||||
//- Orient surface
|
|
||||||
static void orientSurface
|
|
||||||
(
|
|
||||||
triSurface&,
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const Map<labelList>& edgeFacesRest
|
|
||||||
);
|
|
||||||
|
|
||||||
//- Is triangle (given by 3 edges) not fully connected?
|
|
||||||
static bool danglingTriangle
|
|
||||||
(
|
|
||||||
const FixedList<label, 3>& fEdges,
|
|
||||||
const labelList& edgeFace1
|
|
||||||
);
|
|
||||||
|
|
||||||
//- Mark all non-fully connected triangles
|
|
||||||
static label markDanglingTriangles
|
|
||||||
(
|
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
const labelList& edgeFace0,
|
|
||||||
const labelList& edgeFace1,
|
|
||||||
const Map<labelList>& edgeFacesRest,
|
|
||||||
boolList& keepTriangles
|
|
||||||
);
|
|
||||||
|
|
||||||
static triSurface subsetMesh
|
static triSurface subsetMesh
|
||||||
(
|
(
|
||||||
const triSurface& s,
|
const triSurface& s,
|
||||||
@ -392,6 +345,10 @@ class isoSurface
|
|||||||
|
|
||||||
public:
|
public:
|
||||||
|
|
||||||
|
//- Declare friendship with isoSurfaceCell to share some functionality
|
||||||
|
friend class isoSurfaceCell;
|
||||||
|
|
||||||
|
|
||||||
//- Runtime type information
|
//- Runtime type information
|
||||||
TypeName("isoSurface");
|
TypeName("isoSurface");
|
||||||
|
|
||||||
@ -413,18 +370,12 @@ public:
|
|||||||
|
|
||||||
// Member Functions
|
// Member Functions
|
||||||
|
|
||||||
//- For every face original cell in mesh
|
//- For every triangle the original cell in mesh
|
||||||
const labelList& meshCells() const
|
const labelList& meshCells() const
|
||||||
{
|
{
|
||||||
return meshCells_;
|
return meshCells_;
|
||||||
}
|
}
|
||||||
|
|
||||||
//- For every unmerged triangle point the point in the triSurface
|
|
||||||
const labelList& triPointMergeMap() const
|
|
||||||
{
|
|
||||||
return triPointMergeMap_;
|
|
||||||
}
|
|
||||||
|
|
||||||
//- Interpolates cCoords,pCoords. Uses the references to the original
|
//- Interpolates cCoords,pCoords. Uses the references to the original
|
||||||
// fields used to create the iso surface.
|
// fields used to create the iso surface.
|
||||||
template<class Type>
|
template<class Type>
|
||||||
|
|||||||
@ -35,7 +35,7 @@ License
|
|||||||
|
|
||||||
namespace Foam
|
namespace Foam
|
||||||
{
|
{
|
||||||
defineTypeNameAndDebug(isoSurfaceCell, 0);
|
defineTypeNameAndDebug(isoSurfaceCell, 0);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@ -212,8 +212,6 @@ void Foam::isoSurfaceCell::calcCutTypes
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
// Return the two common points between two triangles
|
|
||||||
Foam::labelPair Foam::isoSurfaceCell::findCommonPoints
|
Foam::labelPair Foam::isoSurfaceCell::findCommonPoints
|
||||||
(
|
(
|
||||||
const labelledTri& tri0,
|
const labelledTri& tri0,
|
||||||
@ -250,7 +248,6 @@ Foam::labelPair Foam::isoSurfaceCell::findCommonPoints
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Caculate centre of surface.
|
|
||||||
Foam::point Foam::isoSurfaceCell::calcCentre(const triSurface& s)
|
Foam::point Foam::isoSurfaceCell::calcCentre(const triSurface& s)
|
||||||
{
|
{
|
||||||
vector sum = vector::zero;
|
vector sum = vector::zero;
|
||||||
@ -263,8 +260,6 @@ Foam::point Foam::isoSurfaceCell::calcCentre(const triSurface& s)
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Replace surface (localPoints, localTris) with single point. Returns
|
|
||||||
// point. Destructs arguments.
|
|
||||||
Foam::pointIndexHit Foam::isoSurfaceCell::collapseSurface
|
Foam::pointIndexHit Foam::isoSurfaceCell::collapseSurface
|
||||||
(
|
(
|
||||||
const label cellI,
|
const label cellI,
|
||||||
@ -532,7 +527,6 @@ void Foam::isoSurfaceCell::calcSnappedCc
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Generate triangles for face connected to pointI
|
|
||||||
void Foam::isoSurfaceCell::genPointTris
|
void Foam::isoSurfaceCell::genPointTris
|
||||||
(
|
(
|
||||||
const scalarField& cellValues,
|
const scalarField& cellValues,
|
||||||
@ -605,7 +599,6 @@ void Foam::isoSurfaceCell::genPointTris
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Generate triangle for tet connected to pointI
|
|
||||||
void Foam::isoSurfaceCell::genPointTris
|
void Foam::isoSurfaceCell::genPointTris
|
||||||
(
|
(
|
||||||
const scalarField& pointValues,
|
const scalarField& pointValues,
|
||||||
@ -1053,7 +1046,6 @@ Foam::triSurface Foam::isoSurfaceCell::stitchTriPoints
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Does face use valid vertices?
|
|
||||||
bool Foam::isoSurfaceCell::validTri(const triSurface& surf, const label faceI)
|
bool Foam::isoSurfaceCell::validTri(const triSurface& surf, const label faceI)
|
||||||
{
|
{
|
||||||
// Simple check on indices ok.
|
// Simple check on indices ok.
|
||||||
@ -1222,143 +1214,6 @@ void Foam::isoSurfaceCell::calcAddressing
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
//void Foam::isoSurfaceCell::walkOrientation
|
|
||||||
//(
|
|
||||||
// const triSurface& surf,
|
|
||||||
// const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
// const labelList& edgeFace0,
|
|
||||||
// const labelList& edgeFace1,
|
|
||||||
// const label seedTriI,
|
|
||||||
// labelList& flipState
|
|
||||||
//)
|
|
||||||
//{
|
|
||||||
// // Do walk for consistent orientation.
|
|
||||||
// DynamicList<label> changedFaces(surf.size());
|
|
||||||
//
|
|
||||||
// changedFaces.append(seedTriI);
|
|
||||||
//
|
|
||||||
// while (changedFaces.size())
|
|
||||||
// {
|
|
||||||
// DynamicList<label> newChangedFaces(changedFaces.size());
|
|
||||||
//
|
|
||||||
// forAll(changedFaces, i)
|
|
||||||
// {
|
|
||||||
// label triI = changedFaces[i];
|
|
||||||
// const labelledTri& tri = surf[triI];
|
|
||||||
// const FixedList<label, 3>& fEdges = faceEdges[triI];
|
|
||||||
//
|
|
||||||
// forAll(fEdges, fp)
|
|
||||||
// {
|
|
||||||
// label edgeI = fEdges[fp];
|
|
||||||
//
|
|
||||||
// // my points:
|
|
||||||
// label p0 = tri[fp];
|
|
||||||
// label p1 = tri[tri.fcIndex(fp)];
|
|
||||||
//
|
|
||||||
// label nbrI =
|
|
||||||
// (
|
|
||||||
// edgeFace0[edgeI] != triI
|
|
||||||
// ? edgeFace0[edgeI]
|
|
||||||
// : edgeFace1[edgeI]
|
|
||||||
// );
|
|
||||||
//
|
|
||||||
// if (nbrI != -1 && flipState[nbrI] == -1)
|
|
||||||
// {
|
|
||||||
// const labelledTri& nbrTri = surf[nbrI];
|
|
||||||
//
|
|
||||||
// // nbr points
|
|
||||||
// label nbrFp = findIndex(nbrTri, p0);
|
|
||||||
// label nbrP1 = nbrTri[nbrTri.rcIndex(nbrFp)];
|
|
||||||
//
|
|
||||||
// bool sameOrientation = (p1 == nbrP1);
|
|
||||||
//
|
|
||||||
// if (flipState[triI] == 0)
|
|
||||||
// {
|
|
||||||
// flipState[nbrI] = (sameOrientation ? 0 : 1);
|
|
||||||
// }
|
|
||||||
// else
|
|
||||||
// {
|
|
||||||
// flipState[nbrI] = (sameOrientation ? 1 : 0);
|
|
||||||
// }
|
|
||||||
// newChangedFaces.append(nbrI);
|
|
||||||
// }
|
|
||||||
// }
|
|
||||||
// }
|
|
||||||
//
|
|
||||||
// changedFaces.transfer(newChangedFaces);
|
|
||||||
// }
|
|
||||||
//}
|
|
||||||
//
|
|
||||||
//
|
|
||||||
//void Foam::isoSurfaceCell::orientSurface
|
|
||||||
//(
|
|
||||||
// triSurface& surf,
|
|
||||||
// const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
// const labelList& edgeFace0,
|
|
||||||
// const labelList& edgeFace1,
|
|
||||||
// const Map<labelList>& edgeFacesRest
|
|
||||||
//)
|
|
||||||
//{
|
|
||||||
// // -1 : unvisited
|
|
||||||
// // 0 : leave as is
|
|
||||||
// // 1 : flip
|
|
||||||
// labelList flipState(surf.size(), -1);
|
|
||||||
//
|
|
||||||
// label seedTriI = 0;
|
|
||||||
//
|
|
||||||
// while (true)
|
|
||||||
// {
|
|
||||||
// // Find first unvisited triangle
|
|
||||||
// for
|
|
||||||
// (
|
|
||||||
// ;
|
|
||||||
// seedTriI < surf.size() && flipState[seedTriI] != -1;
|
|
||||||
// seedTriI++
|
|
||||||
// )
|
|
||||||
// {}
|
|
||||||
//
|
|
||||||
// if (seedTriI == surf.size())
|
|
||||||
// {
|
|
||||||
// break;
|
|
||||||
// }
|
|
||||||
//
|
|
||||||
// // Note: Determine orientation of seedTriI?
|
|
||||||
// // for now assume it is ok
|
|
||||||
// flipState[seedTriI] = 0;
|
|
||||||
//
|
|
||||||
// walkOrientation
|
|
||||||
// (
|
|
||||||
// surf,
|
|
||||||
// faceEdges,
|
|
||||||
// edgeFace0,
|
|
||||||
// edgeFace1,
|
|
||||||
// seedTriI,
|
|
||||||
// flipState
|
|
||||||
// );
|
|
||||||
// }
|
|
||||||
//
|
|
||||||
// // Do actual flipping
|
|
||||||
// surf.clearOut();
|
|
||||||
// forAll(surf, triI)
|
|
||||||
// {
|
|
||||||
// if (flipState[triI] == 1)
|
|
||||||
// {
|
|
||||||
// labelledTri tri(surf[triI]);
|
|
||||||
//
|
|
||||||
// surf[triI][0] = tri[0];
|
|
||||||
// surf[triI][1] = tri[2];
|
|
||||||
// surf[triI][2] = tri[1];
|
|
||||||
// }
|
|
||||||
// else if (flipState[triI] == -1)
|
|
||||||
// {
|
|
||||||
// FatalErrorInFunction
|
|
||||||
// << "problem" << abort(FatalError);
|
|
||||||
// }
|
|
||||||
// }
|
|
||||||
//}
|
|
||||||
|
|
||||||
|
|
||||||
// Checks if triangle is connected through edgeI only.
|
|
||||||
bool Foam::isoSurfaceCell::danglingTriangle
|
bool Foam::isoSurfaceCell::danglingTriangle
|
||||||
(
|
(
|
||||||
const FixedList<label, 3>& fEdges,
|
const FixedList<label, 3>& fEdges,
|
||||||
@ -1385,7 +1240,6 @@ bool Foam::isoSurfaceCell::danglingTriangle
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
// Mark triangles to keep. Returns number of dangling triangles.
|
|
||||||
Foam::label Foam::isoSurfaceCell::markDanglingTriangles
|
Foam::label Foam::isoSurfaceCell::markDanglingTriangles
|
||||||
(
|
(
|
||||||
const List<FixedList<label, 3>>& faceEdges,
|
const List<FixedList<label, 3>>& faceEdges,
|
||||||
@ -1605,56 +1459,58 @@ Foam::isoSurfaceCell::isoSurfaceCell
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|
||||||
DynamicList<point> triPoints(nCutCells_);
|
|
||||||
DynamicList<label> triMeshCells(nCutCells_);
|
|
||||||
|
|
||||||
generateTriPoints
|
|
||||||
(
|
|
||||||
cVals,
|
|
||||||
pVals,
|
|
||||||
|
|
||||||
mesh_.cellCentres(),
|
|
||||||
mesh_.points(),
|
|
||||||
|
|
||||||
snappedPoints,
|
|
||||||
snappedCc,
|
|
||||||
snappedPoint,
|
|
||||||
|
|
||||||
triPoints,
|
|
||||||
triMeshCells
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
|
||||||
{
|
{
|
||||||
Pout<< "isoSurfaceCell : generated " << triMeshCells.size()
|
DynamicList<point> triPoints(nCutCells_);
|
||||||
<< " unmerged triangles." << endl;
|
DynamicList<label> triMeshCells(nCutCells_);
|
||||||
}
|
|
||||||
|
|
||||||
// Merge points and compact out non-valid triangles
|
generateTriPoints
|
||||||
labelList triMap; // merged to unmerged triangle
|
|
||||||
triSurface::operator=
|
|
||||||
(
|
|
||||||
stitchTriPoints
|
|
||||||
(
|
(
|
||||||
regularise, // check for duplicate tris
|
cVals,
|
||||||
|
pVals,
|
||||||
|
|
||||||
|
mesh_.cellCentres(),
|
||||||
|
mesh_.points(),
|
||||||
|
|
||||||
|
snappedPoints,
|
||||||
|
snappedCc,
|
||||||
|
snappedPoint,
|
||||||
|
|
||||||
triPoints,
|
triPoints,
|
||||||
triPointMergeMap_, // unmerged to merged point
|
triMeshCells
|
||||||
triMap
|
);
|
||||||
)
|
|
||||||
);
|
|
||||||
|
|
||||||
if (debug)
|
if (debug)
|
||||||
{
|
{
|
||||||
Pout<< "isoSurfaceCell : generated " << triMap.size()
|
Pout<< "isoSurfaceCell : generated " << triMeshCells.size()
|
||||||
<< " merged triangles." << endl;
|
<< " unmerged triangles." << endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
// Merge points and compact out non-valid triangles
|
||||||
|
labelList triMap;
|
||||||
|
triSurface::operator=
|
||||||
|
(
|
||||||
|
stitchTriPoints
|
||||||
|
(
|
||||||
|
regularise, // check for duplicate tris
|
||||||
|
triPoints,
|
||||||
|
triPointMergeMap_, // unmerged to merged point
|
||||||
|
triMap // merged to unmerged triangle
|
||||||
|
)
|
||||||
|
);
|
||||||
|
|
||||||
|
if (debug)
|
||||||
|
{
|
||||||
|
Pout<< "isoSurfaceCell : generated " << triMap.size()
|
||||||
|
<< " merged triangles." << endl;
|
||||||
|
}
|
||||||
|
|
||||||
|
meshCells_.setSize(triMap.size());
|
||||||
|
forAll(triMap, i)
|
||||||
|
{
|
||||||
|
meshCells_[i] = triMeshCells[triMap[i]];
|
||||||
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
meshCells_.setSize(triMap.size());
|
|
||||||
forAll(triMap, i)
|
|
||||||
{
|
|
||||||
meshCells_[i] = triMeshCells[triMap[i]];
|
|
||||||
}
|
|
||||||
|
|
||||||
if (debug)
|
if (debug)
|
||||||
{
|
{
|
||||||
@ -1728,8 +1584,6 @@ Foam::isoSurfaceCell::isoSurfaceCell
|
|||||||
meshCells_ = labelField(meshCells_, subsetTriMap);
|
meshCells_ = labelField(meshCells_, subsetTriMap);
|
||||||
inplaceRenumber(reversePointMap, triPointMergeMap_);
|
inplaceRenumber(reversePointMap, triPointMergeMap_);
|
||||||
}
|
}
|
||||||
|
|
||||||
//orientSurface(*this, faceEdges, edgeFace0, edgeFace1, edgeFacesRest);
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|||||||
@ -139,14 +139,18 @@ class isoSurfaceCell
|
|||||||
const scalarField& pointValues
|
const scalarField& pointValues
|
||||||
);
|
);
|
||||||
|
|
||||||
|
//- Return the two common points between two triangles
|
||||||
static labelPair findCommonPoints
|
static labelPair findCommonPoints
|
||||||
(
|
(
|
||||||
const labelledTri&,
|
const labelledTri&,
|
||||||
const labelledTri&
|
const labelledTri&
|
||||||
);
|
);
|
||||||
|
|
||||||
|
//- Caculate centre of surface.
|
||||||
static point calcCentre(const triSurface&);
|
static point calcCentre(const triSurface&);
|
||||||
|
|
||||||
|
//- Replace surface (localPoints, localTris) with single point.
|
||||||
|
// Returns point. Destroys arguments.
|
||||||
pointIndexHit collapseSurface
|
pointIndexHit collapseSurface
|
||||||
(
|
(
|
||||||
const label cellI,
|
const label cellI,
|
||||||
@ -214,19 +218,15 @@ class isoSurfaceCell
|
|||||||
void generateTriPoints
|
void generateTriPoints
|
||||||
(
|
(
|
||||||
const DynamicList<Type>& snapped,
|
const DynamicList<Type>& snapped,
|
||||||
const scalar isoVal0,
|
|
||||||
const scalar s0,
|
const scalar s0,
|
||||||
const Type& p0,
|
const Type& p0,
|
||||||
const label p0Index,
|
const label p0Index,
|
||||||
const scalar isoVal1,
|
|
||||||
const scalar s1,
|
const scalar s1,
|
||||||
const Type& p1,
|
const Type& p1,
|
||||||
const label p1Index,
|
const label p1Index,
|
||||||
const scalar isoVal2,
|
|
||||||
const scalar s2,
|
const scalar s2,
|
||||||
const Type& p2,
|
const Type& p2,
|
||||||
const label p2Index,
|
const label p2Index,
|
||||||
const scalar isoVal3,
|
|
||||||
const scalar s3,
|
const scalar s3,
|
||||||
const Type& p3,
|
const Type& p3,
|
||||||
const label p3Index,
|
const label p3Index,
|
||||||
@ -271,27 +271,6 @@ class isoSurfaceCell
|
|||||||
Map<labelList>& edgeFacesRest
|
Map<labelList>& edgeFacesRest
|
||||||
) const;
|
) const;
|
||||||
|
|
||||||
////- Determine orientation
|
|
||||||
//static void walkOrientation
|
|
||||||
//(
|
|
||||||
// const triSurface& surf,
|
|
||||||
// const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
// const labelList& edgeFace0,
|
|
||||||
// const labelList& edgeFace1,
|
|
||||||
// const label seedTriI,
|
|
||||||
// labelList& flipState
|
|
||||||
//);
|
|
||||||
|
|
||||||
////- Orient surface
|
|
||||||
//static void orientSurface
|
|
||||||
//(
|
|
||||||
// triSurface&,
|
|
||||||
// const List<FixedList<label, 3>>& faceEdges,
|
|
||||||
// const labelList& edgeFace0,
|
|
||||||
// const labelList& edgeFace1,
|
|
||||||
// const Map<labelList>& edgeFacesRest
|
|
||||||
//);
|
|
||||||
|
|
||||||
//- Is triangle (given by 3 edges) not fully connected?
|
//- Is triangle (given by 3 edges) not fully connected?
|
||||||
static bool danglingTriangle
|
static bool danglingTriangle
|
||||||
(
|
(
|
||||||
@ -317,8 +296,6 @@ class isoSurfaceCell
|
|||||||
labelList& newToOldPoints
|
labelList& newToOldPoints
|
||||||
);
|
);
|
||||||
|
|
||||||
//- Combine all triangles inside a cell into a minimal triangulation
|
|
||||||
void combineCellTriangles();
|
|
||||||
|
|
||||||
public:
|
public:
|
||||||
|
|
||||||
@ -348,24 +325,6 @@ public:
|
|||||||
return meshCells_;
|
return meshCells_;
|
||||||
}
|
}
|
||||||
|
|
||||||
//- For every unmerged triangle point the point in the triSurface
|
|
||||||
const labelList triPointMergeMap() const
|
|
||||||
{
|
|
||||||
return triPointMergeMap_;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
//- Interpolates cCoords,pCoords. Takes the original fields
|
|
||||||
// used to create the iso surface.
|
|
||||||
template<class Type>
|
|
||||||
tmp<Field<Type>> interpolate
|
|
||||||
(
|
|
||||||
const scalarField& cVals,
|
|
||||||
const scalarField& pVals,
|
|
||||||
const Field<Type>& cCoords,
|
|
||||||
const Field<Type>& pCoords
|
|
||||||
) const;
|
|
||||||
|
|
||||||
//- Interpolates cCoords,pCoords.
|
//- Interpolates cCoords,pCoords.
|
||||||
template<class Type>
|
template<class Type>
|
||||||
tmp<Field<Type>> interpolate
|
tmp<Field<Type>> interpolate
|
||||||
|
|||||||
@ -26,6 +26,7 @@ License
|
|||||||
#include "isoSurfaceCell.H"
|
#include "isoSurfaceCell.H"
|
||||||
#include "polyMesh.H"
|
#include "polyMesh.H"
|
||||||
#include "tetMatcher.H"
|
#include "tetMatcher.H"
|
||||||
|
#include "isoSurface.H"
|
||||||
|
|
||||||
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
|
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
|
||||||
|
|
||||||
@ -76,22 +77,18 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
(
|
(
|
||||||
const DynamicList<Type>& snapped,
|
const DynamicList<Type>& snapped,
|
||||||
|
|
||||||
const scalar isoVal0,
|
|
||||||
const scalar s0,
|
const scalar s0,
|
||||||
const Type& p0,
|
const Type& p0,
|
||||||
const label p0Index,
|
const label p0Index,
|
||||||
|
|
||||||
const scalar isoVal1,
|
|
||||||
const scalar s1,
|
const scalar s1,
|
||||||
const Type& p1,
|
const Type& p1,
|
||||||
const label p1Index,
|
const label p1Index,
|
||||||
|
|
||||||
const scalar isoVal2,
|
|
||||||
const scalar s2,
|
const scalar s2,
|
||||||
const Type& p2,
|
const Type& p2,
|
||||||
const label p2Index,
|
const label p2Index,
|
||||||
|
|
||||||
const scalar isoVal3,
|
|
||||||
const scalar s3,
|
const scalar s3,
|
||||||
const Type& p3,
|
const Type& p3,
|
||||||
const label p3Index,
|
const label p3Index,
|
||||||
@ -117,167 +114,203 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
triIndex |= 8;
|
triIndex |= 8;
|
||||||
}
|
}
|
||||||
|
|
||||||
// Form the vertices of the triangles for each case
|
/* Form the vertices of the triangles for each case */
|
||||||
switch (triIndex)
|
switch (triIndex)
|
||||||
{
|
{
|
||||||
case 0x00:
|
case 0x00:
|
||||||
case 0x0F:
|
case 0x0F:
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0E:
|
|
||||||
case 0x01:
|
case 0x01:
|
||||||
|
case 0x0E:
|
||||||
{
|
{
|
||||||
// 0 is common point. Orient such that normal points in positive
|
pts.append
|
||||||
// gradient direction
|
(
|
||||||
if (isoVal0 >= isoVal1)
|
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
if (triIndex == 0x0E)
|
||||||
{
|
{
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index));
|
// Flip normals
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
|
label sz = pts.size();
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index));
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0D:
|
|
||||||
case 0x02:
|
case 0x02:
|
||||||
|
case 0x0D:
|
||||||
{
|
{
|
||||||
// 1 is common point
|
pts.append
|
||||||
if (isoVal1 >= isoVal0)
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
|
||||||
|
if (triIndex == 0x0D)
|
||||||
{
|
{
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index));
|
// Flip normals
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
|
label sz = pts.size();
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s0,p0,p0Index));
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0C:
|
|
||||||
case 0x03:
|
case 0x03:
|
||||||
|
case 0x0C:
|
||||||
{
|
{
|
||||||
Type s02 = generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
|
Type p0p2 =
|
||||||
Type s13 = generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
|
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index);
|
||||||
|
Type p1p3 =
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index);
|
||||||
|
|
||||||
if (isoVal0 >= isoVal3)
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
pts.append(p1p3);
|
||||||
|
pts.append(p0p2);
|
||||||
|
|
||||||
|
pts.append(p1p3);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
pts.append(p0p2);
|
||||||
|
|
||||||
|
if (triIndex == 0x0C)
|
||||||
{
|
{
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
// Flip normals
|
||||||
pts.append(s02);
|
label sz = pts.size();
|
||||||
pts.append(s13);
|
Swap(pts[sz-5], pts[sz-4]);
|
||||||
pts.append(s13);
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
|
||||||
pts.append(s02);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(s02);
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
|
||||||
pts.append(s13);
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
|
||||||
pts.append(s13);
|
|
||||||
pts.append(s02);
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0B:
|
|
||||||
case 0x04:
|
case 0x04:
|
||||||
|
case 0x0B:
|
||||||
{
|
{
|
||||||
// 2 is common point
|
pts.append
|
||||||
if (isoVal2 >= isoVal0)
|
(
|
||||||
|
generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
|
||||||
|
if (triIndex == 0x0B)
|
||||||
{
|
{
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index));
|
// Flip normals
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index));
|
label sz = pts.size();
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index));
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s1,p1,p1Index));
|
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s0,p0,p0Index));
|
|
||||||
pts.append(generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index));
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0A:
|
|
||||||
case 0x05:
|
case 0x05:
|
||||||
|
case 0x0A:
|
||||||
{
|
{
|
||||||
Type s01 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
|
Type p0p1 =
|
||||||
Type s23 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
|
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
|
||||||
|
Type p2p3 =
|
||||||
|
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
|
||||||
|
|
||||||
if (isoVal3 >= isoVal0)
|
pts.append(p0p1);
|
||||||
|
pts.append(p2p3);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
|
||||||
|
pts.append(p0p1);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
pts.append(p2p3);
|
||||||
|
|
||||||
|
if (triIndex == 0x0A)
|
||||||
{
|
{
|
||||||
pts.append(s01);
|
// Flip normals
|
||||||
pts.append(s23);
|
label sz = pts.size();
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
Swap(pts[sz-5], pts[sz-4]);
|
||||||
pts.append(s01);
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
|
||||||
pts.append(s23);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(s23);
|
|
||||||
pts.append(s01);
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s3,p3,p3Index));
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s2,p2,p2Index));
|
|
||||||
pts.append(s01);
|
|
||||||
pts.append(s23);
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x09:
|
|
||||||
case 0x06:
|
case 0x06:
|
||||||
|
case 0x09:
|
||||||
{
|
{
|
||||||
Type s01 = generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
|
Type p0p1 =
|
||||||
Type s23 = generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
|
generatePoint(snapped,s0,p0,p0Index,s1,p1,p1Index);
|
||||||
|
Type p2p3 =
|
||||||
|
generatePoint(snapped,s2,p2,p2Index,s3,p3,p3Index);
|
||||||
|
|
||||||
if (isoVal3 >= isoVal1)
|
pts.append(p0p1);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index)
|
||||||
|
);
|
||||||
|
pts.append(p2p3);
|
||||||
|
|
||||||
|
pts.append(p0p1);
|
||||||
|
pts.append(p2p3);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
|
||||||
|
if (triIndex == 0x09)
|
||||||
{
|
{
|
||||||
pts.append(s01);
|
// Flip normals
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
|
label sz = pts.size();
|
||||||
pts.append(s23);
|
Swap(pts[sz-5], pts[sz-4]);
|
||||||
pts.append(s01);
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
|
|
||||||
pts.append(s23);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(generatePoint(snapped,s1,p1,p1Index,s3,p3,p3Index));
|
|
||||||
pts.append(s01);
|
|
||||||
pts.append(s23);
|
|
||||||
pts.append(generatePoint(snapped,s0,p0,p0Index,s2,p2,p2Index));
|
|
||||||
pts.append(s01);
|
|
||||||
pts.append(s23);
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x07:
|
|
||||||
case 0x08:
|
case 0x08:
|
||||||
|
case 0x07:
|
||||||
{
|
{
|
||||||
// 3 is common point
|
pts.append
|
||||||
if (isoVal3 >= isoVal0)
|
(
|
||||||
|
generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index)
|
||||||
|
);
|
||||||
|
pts.append
|
||||||
|
(
|
||||||
|
generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index)
|
||||||
|
);
|
||||||
|
if (triIndex == 0x07)
|
||||||
{
|
{
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index));
|
// Flip normals
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index));
|
label sz = pts.size();
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index));
|
Swap(pts[sz-2], pts[sz-1]);
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s2,p2,p2Index));
|
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s0,p0,p0Index));
|
|
||||||
pts.append(generatePoint(snapped,s3,p3,p3Index,s1,p1,p1Index));
|
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
break;
|
break;
|
||||||
@ -341,22 +374,18 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
(
|
(
|
||||||
snappedPoints,
|
snappedPoints,
|
||||||
|
|
||||||
pVals_[f0[1]],
|
|
||||||
pVals[f0[1]],
|
pVals[f0[1]],
|
||||||
pCoords[f0[1]],
|
pCoords[f0[1]],
|
||||||
snappedPoint[f0[1]],
|
snappedPoint[f0[1]],
|
||||||
|
|
||||||
pVals_[f0[0]],
|
|
||||||
pVals[f0[0]],
|
pVals[f0[0]],
|
||||||
pCoords[f0[0]],
|
pCoords[f0[0]],
|
||||||
snappedPoint[f0[0]],
|
snappedPoint[f0[0]],
|
||||||
|
|
||||||
pVals_[f0[2]],
|
|
||||||
pVals[f0[2]],
|
pVals[f0[2]],
|
||||||
pCoords[f0[2]],
|
pCoords[f0[2]],
|
||||||
snappedPoint[f0[2]],
|
snappedPoint[f0[2]],
|
||||||
|
|
||||||
pVals_[oppositeI],
|
|
||||||
pVals[oppositeI],
|
pVals[oppositeI],
|
||||||
pCoords[oppositeI],
|
pCoords[oppositeI],
|
||||||
snappedPoint[oppositeI],
|
snappedPoint[oppositeI],
|
||||||
@ -370,22 +399,18 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
(
|
(
|
||||||
snappedPoints,
|
snappedPoints,
|
||||||
|
|
||||||
pVals_[f0[0]],
|
|
||||||
pVals[f0[0]],
|
pVals[f0[0]],
|
||||||
pCoords[f0[0]],
|
pCoords[f0[0]],
|
||||||
snappedPoint[f0[0]],
|
snappedPoint[f0[0]],
|
||||||
|
|
||||||
pVals_[f0[1]],
|
|
||||||
pVals[f0[1]],
|
pVals[f0[1]],
|
||||||
pCoords[f0[1]],
|
pCoords[f0[1]],
|
||||||
snappedPoint[f0[1]],
|
snappedPoint[f0[1]],
|
||||||
|
|
||||||
pVals_[f0[2]],
|
|
||||||
pVals[f0[2]],
|
pVals[f0[2]],
|
||||||
pCoords[f0[2]],
|
pCoords[f0[2]],
|
||||||
snappedPoint[f0[2]],
|
snappedPoint[f0[2]],
|
||||||
|
|
||||||
pVals_[oppositeI],
|
|
||||||
pVals[oppositeI],
|
pVals[oppositeI],
|
||||||
pCoords[oppositeI],
|
pCoords[oppositeI],
|
||||||
snappedPoint[oppositeI],
|
snappedPoint[oppositeI],
|
||||||
@ -411,7 +436,6 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
}
|
}
|
||||||
|
|
||||||
label fp = f.fcIndex(fp0);
|
label fp = f.fcIndex(fp0);
|
||||||
|
|
||||||
for (label i = 2; i < f.size(); i++)
|
for (label i = 2; i < f.size(); i++)
|
||||||
{
|
{
|
||||||
label nextFp = f.fcIndex(fp);
|
label nextFp = f.fcIndex(fp);
|
||||||
@ -425,22 +449,18 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
(
|
(
|
||||||
snappedPoints,
|
snappedPoints,
|
||||||
|
|
||||||
pVals_[tri[1]],
|
|
||||||
pVals[tri[1]],
|
pVals[tri[1]],
|
||||||
pCoords[tri[1]],
|
pCoords[tri[1]],
|
||||||
snappedPoint[tri[1]],
|
snappedPoint[tri[1]],
|
||||||
|
|
||||||
pVals_[tri[0]],
|
|
||||||
pVals[tri[0]],
|
pVals[tri[0]],
|
||||||
pCoords[tri[0]],
|
pCoords[tri[0]],
|
||||||
snappedPoint[tri[0]],
|
snappedPoint[tri[0]],
|
||||||
|
|
||||||
pVals_[tri[2]],
|
|
||||||
pVals[tri[2]],
|
pVals[tri[2]],
|
||||||
pCoords[tri[2]],
|
pCoords[tri[2]],
|
||||||
snappedPoint[tri[2]],
|
snappedPoint[tri[2]],
|
||||||
|
|
||||||
cVals_[cellI],
|
|
||||||
cVals[cellI],
|
cVals[cellI],
|
||||||
cCoords[cellI],
|
cCoords[cellI],
|
||||||
snappedCc[cellI],
|
snappedCc[cellI],
|
||||||
@ -454,22 +474,18 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
(
|
(
|
||||||
snappedPoints,
|
snappedPoints,
|
||||||
|
|
||||||
pVals_[tri[0]],
|
|
||||||
pVals[tri[0]],
|
pVals[tri[0]],
|
||||||
pCoords[tri[0]],
|
pCoords[tri[0]],
|
||||||
snappedPoint[tri[0]],
|
snappedPoint[tri[0]],
|
||||||
|
|
||||||
pVals_[tri[1]],
|
|
||||||
pVals[tri[1]],
|
pVals[tri[1]],
|
||||||
pCoords[tri[1]],
|
pCoords[tri[1]],
|
||||||
snappedPoint[tri[1]],
|
snappedPoint[tri[1]],
|
||||||
|
|
||||||
pVals_[tri[2]],
|
|
||||||
pVals[tri[2]],
|
pVals[tri[2]],
|
||||||
pCoords[tri[2]],
|
pCoords[tri[2]],
|
||||||
snappedPoint[tri[2]],
|
snappedPoint[tri[2]],
|
||||||
|
|
||||||
cVals_[cellI],
|
|
||||||
cVals[cellI],
|
cVals[cellI],
|
||||||
cCoords[cellI],
|
cCoords[cellI],
|
||||||
snappedCc[cellI],
|
snappedCc[cellI],
|
||||||
@ -495,7 +511,7 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
|
|
||||||
if (countNotFoundTets > 0)
|
if (countNotFoundTets > 0)
|
||||||
{
|
{
|
||||||
WarningInFunction
|
WarningIn("Foam::isoSurfaceCell::generateTriPoints(..)")
|
||||||
<< "Could not find " << countNotFoundTets
|
<< "Could not find " << countNotFoundTets
|
||||||
<< " tet base points, which may lead to inverted triangles."
|
<< " tet base points, which may lead to inverted triangles."
|
||||||
<< endl;
|
<< endl;
|
||||||
@ -507,68 +523,14 @@ void Foam::isoSurfaceCell::generateTriPoints
|
|||||||
|
|
||||||
|
|
||||||
template<class Type>
|
template<class Type>
|
||||||
Foam::tmp<Foam::Field<Type>>
|
Foam::tmp<Foam::Field<Type> >
|
||||||
Foam::isoSurfaceCell::interpolate
|
|
||||||
(
|
|
||||||
const scalarField& cVals,
|
|
||||||
const scalarField& pVals,
|
|
||||||
const Field<Type>& cCoords,
|
|
||||||
const Field<Type>& pCoords
|
|
||||||
) const
|
|
||||||
{
|
|
||||||
DynamicList<Type> triPoints(nCutCells_);
|
|
||||||
DynamicList<label> triMeshCells(nCutCells_);
|
|
||||||
|
|
||||||
// Dummy snap data
|
|
||||||
DynamicList<Type> snappedPoints;
|
|
||||||
labelList snappedCc(mesh_.nCells(), -1);
|
|
||||||
labelList snappedPoint(mesh_.nPoints(), -1);
|
|
||||||
|
|
||||||
|
|
||||||
generateTriPoints
|
|
||||||
(
|
|
||||||
cVals,
|
|
||||||
pVals,
|
|
||||||
|
|
||||||
cCoords,
|
|
||||||
pCoords,
|
|
||||||
|
|
||||||
snappedPoints,
|
|
||||||
snappedCc,
|
|
||||||
snappedPoint,
|
|
||||||
|
|
||||||
triPoints,
|
|
||||||
triMeshCells
|
|
||||||
);
|
|
||||||
|
|
||||||
|
|
||||||
// One value per point
|
|
||||||
tmp<Field<Type>> tvalues(new Field<Type>(points().size()));
|
|
||||||
Field<Type>& values = tvalues();
|
|
||||||
|
|
||||||
forAll(triPoints, i)
|
|
||||||
{
|
|
||||||
label mergedPointI = triPointMergeMap_[i];
|
|
||||||
|
|
||||||
if (mergedPointI >= 0)
|
|
||||||
{
|
|
||||||
values[mergedPointI] = triPoints[i];
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
return tvalues;
|
|
||||||
}
|
|
||||||
|
|
||||||
|
|
||||||
template<class Type>
|
|
||||||
Foam::tmp<Foam::Field<Type>>
|
|
||||||
Foam::isoSurfaceCell::interpolate
|
Foam::isoSurfaceCell::interpolate
|
||||||
(
|
(
|
||||||
const Field<Type>& cCoords,
|
const Field<Type>& cCoords,
|
||||||
const Field<Type>& pCoords
|
const Field<Type>& pCoords
|
||||||
) const
|
) const
|
||||||
{
|
{
|
||||||
DynamicList<Type> triPoints(nCutCells_);
|
DynamicList<Type> triPoints(3*nCutCells_);
|
||||||
DynamicList<label> triMeshCells(nCutCells_);
|
DynamicList<label> triMeshCells(nCutCells_);
|
||||||
|
|
||||||
// Dummy snap data
|
// Dummy snap data
|
||||||
@ -576,7 +538,6 @@ Foam::isoSurfaceCell::interpolate
|
|||||||
labelList snappedCc(mesh_.nCells(), -1);
|
labelList snappedCc(mesh_.nCells(), -1);
|
||||||
labelList snappedPoint(mesh_.nPoints(), -1);
|
labelList snappedPoint(mesh_.nPoints(), -1);
|
||||||
|
|
||||||
|
|
||||||
generateTriPoints
|
generateTriPoints
|
||||||
(
|
(
|
||||||
cVals_,
|
cVals_,
|
||||||
@ -593,22 +554,12 @@ Foam::isoSurfaceCell::interpolate
|
|||||||
triMeshCells
|
triMeshCells
|
||||||
);
|
);
|
||||||
|
|
||||||
|
return isoSurface::interpolate
|
||||||
// One value per point
|
(
|
||||||
tmp<Field<Type>> tvalues(new Field<Type>(points().size()));
|
points().size(),
|
||||||
Field<Type>& values = tvalues.ref();
|
triPointMergeMap_,
|
||||||
|
triPoints
|
||||||
forAll(triPoints, i)
|
);
|
||||||
{
|
|
||||||
label mergedPointI = triPointMergeMap_[i];
|
|
||||||
|
|
||||||
if (mergedPointI >= 0)
|
|
||||||
{
|
|
||||||
values[mergedPointI] = triPoints[i];
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
return tvalues;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
@ -120,7 +120,7 @@ Foam::isoSurface::adaptPatchFields
|
|||||||
{
|
{
|
||||||
fvPatchField<Type>& pfld = const_cast<fvPatchField<Type>&>
|
fvPatchField<Type>& pfld = const_cast<fvPatchField<Type>&>
|
||||||
(
|
(
|
||||||
fld.boundaryField()[patchI]
|
sliceFld.boundaryField()[patchI]
|
||||||
);
|
);
|
||||||
|
|
||||||
const scalarField& w = mesh.weights().boundaryField()[patchI];
|
const scalarField& w = mesh.weights().boundaryField()[patchI];
|
||||||
@ -240,8 +240,8 @@ void Foam::isoSurface::generateTriPoints
|
|||||||
case 0x0F:
|
case 0x0F:
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0E:
|
|
||||||
case 0x01:
|
case 0x01:
|
||||||
|
case 0x0E:
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1)
|
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1)
|
||||||
@ -254,10 +254,16 @@ void Foam::isoSurface::generateTriPoints
|
|||||||
(
|
(
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
||||||
);
|
);
|
||||||
|
if (triIndex == 0x0E)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0D:
|
|
||||||
case 0x02:
|
case 0x02:
|
||||||
|
case 0x0D:
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s0,p0,hasSnap0,snapP0)
|
generatePoint(s1,p1,hasSnap1,snapP1,s0,p0,hasSnap0,snapP0)
|
||||||
@ -270,97 +276,133 @@ void Foam::isoSurface::generateTriPoints
|
|||||||
(
|
(
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
||||||
);
|
);
|
||||||
|
if (triIndex == 0x0D)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0C:
|
|
||||||
case 0x03:
|
case 0x03:
|
||||||
{
|
case 0x0C:
|
||||||
Type tp1 =
|
{
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2);
|
Type p0p2 =
|
||||||
Type tp2 =
|
generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2);
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3);
|
Type p1p3 =
|
||||||
|
generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3);
|
||||||
|
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
||||||
);
|
);
|
||||||
points.append(tp1);
|
points.append(p1p3);
|
||||||
points.append(tp2);
|
points.append(p0p2);
|
||||||
points.append(tp2);
|
|
||||||
points.append
|
points.append(p1p3);
|
||||||
(
|
points.append
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
(
|
||||||
);
|
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
||||||
points.append(tp1);
|
);
|
||||||
}
|
points.append(p0p2);
|
||||||
|
}
|
||||||
|
if (triIndex == 0x0C)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-5], points[sz-4]);
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0B:
|
|
||||||
case 0x04:
|
case 0x04:
|
||||||
{
|
case 0x0B:
|
||||||
points.append
|
{
|
||||||
(
|
points.append
|
||||||
generatePoint(s2,p2,hasSnap2,snapP2,s0,p0,hasSnap0,snapP0)
|
(
|
||||||
);
|
generatePoint(s2,p2,hasSnap2,snapP2,s0,p0,hasSnap0,snapP0)
|
||||||
points.append
|
);
|
||||||
(
|
points.append
|
||||||
generatePoint(s2,p2,hasSnap2,snapP2,s1,p1,hasSnap1,snapP1)
|
(
|
||||||
);
|
generatePoint(s2,p2,hasSnap2,snapP2,s1,p1,hasSnap1,snapP1)
|
||||||
points.append
|
);
|
||||||
(
|
points.append
|
||||||
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3)
|
(
|
||||||
);
|
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3)
|
||||||
}
|
);
|
||||||
|
}
|
||||||
|
if (triIndex == 0x0B)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x0A:
|
|
||||||
case 0x05:
|
case 0x05:
|
||||||
{
|
case 0x0A:
|
||||||
Type tp0 =
|
{
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
|
Type p0p1 =
|
||||||
Type tp1 =
|
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
|
||||||
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
|
Type p2p3 =
|
||||||
|
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
|
||||||
|
|
||||||
points.append(tp0);
|
points.append(p0p1);
|
||||||
points.append(tp1);
|
points.append(p2p3);
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
generatePoint(s0,p0,hasSnap0,snapP0,s3,p3,hasSnap3,snapP3)
|
||||||
);
|
);
|
||||||
points.append(tp0);
|
|
||||||
points.append
|
points.append(p0p1);
|
||||||
(
|
points.append
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
(
|
||||||
);
|
generatePoint(s1,p1,hasSnap1,snapP1,s2,p2,hasSnap2,snapP2)
|
||||||
points.append(tp1);
|
);
|
||||||
}
|
points.append(p2p3);
|
||||||
|
}
|
||||||
|
if (triIndex == 0x0A)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-5], points[sz-4]);
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x09:
|
|
||||||
case 0x06:
|
case 0x06:
|
||||||
{
|
case 0x09:
|
||||||
Type tp0 =
|
{
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
|
Type p0p1 =
|
||||||
Type tp1 =
|
generatePoint(s0,p0,hasSnap0,snapP0,s1,p1,hasSnap1,snapP1);
|
||||||
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
|
Type p2p3 =
|
||||||
|
generatePoint(s2,p2,hasSnap2,snapP2,s3,p3,hasSnap3,snapP3);
|
||||||
|
|
||||||
points.append(tp0);
|
points.append(p0p1);
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
|
generatePoint(s1,p1,hasSnap1,snapP1,s3,p3,hasSnap3,snapP3)
|
||||||
);
|
);
|
||||||
points.append(tp1);
|
points.append(p2p3);
|
||||||
points.append(tp0);
|
|
||||||
points.append
|
points.append(p0p1);
|
||||||
(
|
points.append(p2p3);
|
||||||
generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
|
points.append
|
||||||
);
|
(
|
||||||
points.append(tp1);
|
generatePoint(s0,p0,hasSnap0,snapP0,s2,p2,hasSnap2,snapP2)
|
||||||
}
|
);
|
||||||
|
}
|
||||||
|
if (triIndex == 0x09)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-5], points[sz-4]);
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
|
|
||||||
case 0x07:
|
|
||||||
case 0x08:
|
case 0x08:
|
||||||
|
case 0x07:
|
||||||
points.append
|
points.append
|
||||||
(
|
(
|
||||||
generatePoint(s3,p3,hasSnap3,snapP3,s0,p0,hasSnap0,snapP0)
|
generatePoint(s3,p3,hasSnap3,snapP3,s0,p0,hasSnap0,snapP0)
|
||||||
@ -373,6 +415,12 @@ void Foam::isoSurface::generateTriPoints
|
|||||||
(
|
(
|
||||||
generatePoint(s3,p3,hasSnap3,snapP3,s1,p1,hasSnap1,snapP1)
|
generatePoint(s3,p3,hasSnap3,snapP3,s1,p1,hasSnap1,snapP1)
|
||||||
);
|
);
|
||||||
|
if (triIndex == 0x07)
|
||||||
|
{
|
||||||
|
// Flip normals
|
||||||
|
label sz = points.size();
|
||||||
|
Swap(points[sz-2], points[sz-1]);
|
||||||
|
}
|
||||||
break;
|
break;
|
||||||
}
|
}
|
||||||
}
|
}
|
||||||
@ -681,16 +729,45 @@ void Foam::isoSurface::generateTriPoints
|
|||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
//template<class Type>
|
template<class Type>
|
||||||
//Foam::tmp<Foam::Field<Type>>
|
Foam::tmp<Foam::Field<Type>>
|
||||||
//Foam::isoSurface::sample(const Field<Type>& vField) const
|
Foam::isoSurface::interpolate
|
||||||
//{
|
(
|
||||||
// return tmp<Field<Type>>(new Field<Type>(vField, meshCells()));
|
const label nPoints,
|
||||||
//}
|
const labelList& triPointMergeMap,
|
||||||
|
const DynamicList<Type>& unmergedValues
|
||||||
|
)
|
||||||
|
{
|
||||||
|
// One value per point
|
||||||
|
tmp<Field<Type> > tvalues(new Field<Type>(nPoints, pTraits<Type>::zero));
|
||||||
|
Field<Type>& values = tvalues.ref();
|
||||||
|
labelList nValues(values.size(), 0);
|
||||||
|
|
||||||
|
forAll(unmergedValues, i)
|
||||||
|
{
|
||||||
|
label mergedPointI = triPointMergeMap[i];
|
||||||
|
|
||||||
|
if (mergedPointI >= 0)
|
||||||
|
{
|
||||||
|
values[mergedPointI] += unmergedValues[i];
|
||||||
|
nValues[mergedPointI]++;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
forAll(values, i)
|
||||||
|
{
|
||||||
|
if (nValues[i] > 0)
|
||||||
|
{
|
||||||
|
values[i] /= scalar(nValues[i]);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
return tvalues;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
template<class Type>
|
template<class Type>
|
||||||
Foam::tmp<Foam::Field<Type>>
|
Foam::tmp<Foam::Field<Type> >
|
||||||
Foam::isoSurface::interpolate
|
Foam::isoSurface::interpolate
|
||||||
(
|
(
|
||||||
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
|
const GeometricField<Type, fvPatchField, volMesh>& cCoords,
|
||||||
@ -707,7 +784,7 @@ Foam::isoSurface::interpolate
|
|||||||
>> c2(adaptPatchFields(cCoords));
|
>> c2(adaptPatchFields(cCoords));
|
||||||
|
|
||||||
|
|
||||||
DynamicList<Type> triPoints(nCutCells_);
|
DynamicList<Type> triPoints(3*nCutCells_);
|
||||||
DynamicList<label> triMeshCells(nCutCells_);
|
DynamicList<label> triMeshCells(nCutCells_);
|
||||||
|
|
||||||
// Dummy snap data
|
// Dummy snap data
|
||||||
@ -731,52 +808,12 @@ Foam::isoSurface::interpolate
|
|||||||
triMeshCells
|
triMeshCells
|
||||||
);
|
);
|
||||||
|
|
||||||
|
return interpolate
|
||||||
// One value per point
|
|
||||||
tmp<Field<Type>> tvalues
|
|
||||||
(
|
(
|
||||||
new Field<Type>(points().size(), pTraits<Type>::zero)
|
points().size(),
|
||||||
|
triPointMergeMap_,
|
||||||
|
triPoints
|
||||||
);
|
);
|
||||||
Field<Type>& values = tvalues.ref();
|
|
||||||
labelList nValues(values.size(), 0);
|
|
||||||
|
|
||||||
forAll(triPoints, i)
|
|
||||||
{
|
|
||||||
label mergedPointI = triPointMergeMap_[i];
|
|
||||||
|
|
||||||
if (mergedPointI >= 0)
|
|
||||||
{
|
|
||||||
values[mergedPointI] += triPoints[i];
|
|
||||||
nValues[mergedPointI]++;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
if (debug)
|
|
||||||
{
|
|
||||||
Pout<< "nValues:" << values.size() << endl;
|
|
||||||
label nMult = 0;
|
|
||||||
forAll(nValues, i)
|
|
||||||
{
|
|
||||||
if (nValues[i] == 0)
|
|
||||||
{
|
|
||||||
FatalErrorInFunction
|
|
||||||
<< "point:" << i << " nValues:" << nValues[i]
|
|
||||||
<< abort(FatalError);
|
|
||||||
}
|
|
||||||
else if (nValues[i] > 1)
|
|
||||||
{
|
|
||||||
nMult++;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
Pout<< "Of which mult:" << nMult << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
forAll(values, i)
|
|
||||||
{
|
|
||||||
values[i] /= scalar(nValues[i]);
|
|
||||||
}
|
|
||||||
|
|
||||||
return tvalues;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
|
|||||||
@ -124,18 +124,52 @@ void Foam::sampledIsoSurface::getIsoFields() const
|
|||||||
// Get pointField
|
// Get pointField
|
||||||
// ~~~~~~~~~~~~~~
|
// ~~~~~~~~~~~~~~
|
||||||
|
|
||||||
|
// In case of multiple iso values we don't want to calculate multiple e.g.
|
||||||
|
// "volPointInterpolate(p)" so register it and re-use it. This is the
|
||||||
|
// same as the 'cache' functionality from volPointInterpolate but
|
||||||
|
// unfortunately that one does not guarantee that the field pointer
|
||||||
|
// remain: e.g. some other functionObject might delete the cached version.
|
||||||
|
// (volPointInterpolation::interpolate with cache=false deletes any
|
||||||
|
// registered one or if mesh.changing())
|
||||||
|
|
||||||
if (!subMeshPtr_.valid())
|
if (!subMeshPtr_.valid())
|
||||||
{
|
{
|
||||||
word pointFldName = "volPointInterpolate(" + isoField_ + ')';
|
const word pointFldName =
|
||||||
|
"volPointInterpolate_"
|
||||||
|
+ type()
|
||||||
|
+ "("
|
||||||
|
+ isoField_
|
||||||
|
+ ')';
|
||||||
|
|
||||||
if (fvm.foundObject<pointScalarField>(pointFldName))
|
if (fvm.foundObject<pointScalarField>(pointFldName))
|
||||||
{
|
{
|
||||||
if (debug)
|
if (debug)
|
||||||
{
|
{
|
||||||
InfoInFunction
|
InfoInFunction
|
||||||
<< "Lookup pointField " << pointFldName << endl;
|
<< "lookup pointField " << pointFldName << endl;
|
||||||
}
|
}
|
||||||
pointFieldPtr_ = &fvm.lookupObject<pointScalarField>(pointFldName);
|
const pointScalarField& pfld = fvm.lookupObject<pointScalarField>
|
||||||
|
(
|
||||||
|
pointFldName
|
||||||
|
);
|
||||||
|
|
||||||
|
if (!pfld.upToDate(*volFieldPtr_))
|
||||||
|
{
|
||||||
|
if (debug)
|
||||||
|
{
|
||||||
|
InfoInFunction
|
||||||
|
<< "updating pointField "
|
||||||
|
<< pointFldName << endl;
|
||||||
|
}
|
||||||
|
// Update the interpolated value
|
||||||
|
volPointInterpolation::New(fvm).interpolate
|
||||||
|
(
|
||||||
|
*volFieldPtr_,
|
||||||
|
const_cast<pointScalarField&>(pfld)
|
||||||
|
);
|
||||||
|
}
|
||||||
|
|
||||||
|
pointFieldPtr_ = &pfld;
|
||||||
}
|
}
|
||||||
else
|
else
|
||||||
{
|
{
|
||||||
@ -149,27 +183,20 @@ void Foam::sampledIsoSurface::getIsoFields() const
|
|||||||
<< endl;
|
<< endl;
|
||||||
}
|
}
|
||||||
|
|
||||||
if
|
// Interpolate without cache. Note that we're registering it
|
||||||
|
// below so next time round it goes into the condition
|
||||||
|
// above.
|
||||||
|
tmp<pointScalarField> tpfld
|
||||||
(
|
(
|
||||||
storedPointFieldPtr_.empty()
|
volPointInterpolation::New(fvm).interpolate
|
||||||
|| (fvm.time().timeName() != storedPointFieldPtr_().instance())
|
|
||||||
)
|
|
||||||
{
|
|
||||||
if (debug)
|
|
||||||
{
|
|
||||||
InfoInFunction
|
|
||||||
<< "Interpolating volField " << volFieldPtr_->name()
|
|
||||||
<< " to get pointField " << pointFldName << endl;
|
|
||||||
}
|
|
||||||
|
|
||||||
storedPointFieldPtr_.reset
|
|
||||||
(
|
(
|
||||||
volPointInterpolation::New(fvm)
|
*volFieldPtr_,
|
||||||
.interpolate(*volFieldPtr_).ptr()
|
pointFldName,
|
||||||
);
|
false
|
||||||
storedPointFieldPtr_->checkOut();
|
)
|
||||||
pointFieldPtr_ = storedPointFieldPtr_.operator->();
|
);
|
||||||
}
|
pointFieldPtr_ = tpfld.ptr();
|
||||||
|
const_cast<pointScalarField*>(pointFieldPtr_)->store();
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@ -427,7 +454,6 @@ Foam::sampledIsoSurface::sampledIsoSurface
|
|||||||
prevTimeIndex_(-1),
|
prevTimeIndex_(-1),
|
||||||
storedVolFieldPtr_(NULL),
|
storedVolFieldPtr_(NULL),
|
||||||
volFieldPtr_(NULL),
|
volFieldPtr_(NULL),
|
||||||
storedPointFieldPtr_(NULL),
|
|
||||||
pointFieldPtr_(NULL)
|
pointFieldPtr_(NULL)
|
||||||
{
|
{
|
||||||
if (!sampledSurface::interpolate())
|
if (!sampledSurface::interpolate())
|
||||||
|
|||||||
@ -94,7 +94,6 @@ class sampledIsoSurface
|
|||||||
mutable const volScalarField* volFieldPtr_;
|
mutable const volScalarField* volFieldPtr_;
|
||||||
|
|
||||||
//- Cached pointfield
|
//- Cached pointfield
|
||||||
mutable autoPtr<pointScalarField> storedPointFieldPtr_;
|
|
||||||
mutable const pointScalarField* pointFieldPtr_;
|
mutable const pointScalarField* pointFieldPtr_;
|
||||||
|
|
||||||
// And on subsetted mesh
|
// And on subsetted mesh
|
||||||
|
|||||||
Reference in New Issue
Block a user