songtianlun/openfoam
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

Browse Source
This commit is contained in:
mattijs 2013-09-20 17:36:11 +01:00
commit 2a3c54ba8f
19 changed files with 1316 additions and 1590 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

65
applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
View File

@ -95,24 +95,56 @@ int main(int argc, char *argv[])
const bool collapseFaces = args.optionFound("collapseFaces");
const bool collapseFaceZone = args.optionFound("collapseFaceZone");
if (collapseFaces && collapseFaceZone)
{
FatalErrorIn("main(int, char*[])")
<< "Both face zone collapsing and face collapsing have been"
<< "selected. Choose only one of:" << nl
<< " -collapseFaces" << nl
<< " -collapseFaceZone <faceZoneName>"
<< abort(FatalError);
}
labelIOList pointPriority
(
IOobject
(
"pointPriority",
runTime.timeName(),
runTime,
IOobject::READ_IF_PRESENT,
IOobject::AUTO_WRITE
),
labelList(mesh.nPoints(), labelMin)
);
forAll(timeDirs, timeI)
{
runTime.setTime(timeDirs[timeI], timeI);
Info<< "Time = " << runTime.timeName() << endl;
polyMeshFilter meshFilter(mesh);
autoPtr<polyMeshFilter> meshFilterPtr;
// newMesh will be empty until it is filtered
const autoPtr<fvMesh>& newMesh = meshFilter.filteredMesh();
label nBadFaces = 0;
// Filter small edges only. This reduces the number of faces so that
// the face filtering is sped up.
label nBadFaces = meshFilter.filterEdges(0);
{
polyTopoChange meshMod(newMesh);
meshFilterPtr.set(new polyMeshFilter(mesh, pointPriority));
polyMeshFilter& meshFilter = meshFilterPtr();
meshMod.changeMesh(mesh, false);
// newMesh will be empty until it is filtered
const autoPtr<fvMesh>& newMesh = meshFilter.filteredMesh();
// Filter small edges only. This reduces the number of faces so that
// the face filtering is sped up.
nBadFaces = meshFilter.filterEdges(0);
{
polyTopoChange meshMod(newMesh);
meshMod.changeMesh(mesh, false);
}
pointPriority = meshFilter.pointPriority();
}
if (collapseFaceZone)
@ -121,18 +153,30 @@ int main(int argc, char *argv[])
const faceZone& fZone = mesh.faceZones()[faceZoneName];
meshFilterPtr.reset(new polyMeshFilter(mesh, pointPriority));
polyMeshFilter& meshFilter = meshFilterPtr();
const autoPtr<fvMesh>& newMesh = meshFilter.filteredMesh();
// Filter faces. Pass in the number of bad faces that are present
// from the previous edge filtering to use as a stopping criterion.
meshFilter.filterFaceZone(fZone);
meshFilter.filter(fZone);
{
polyTopoChange meshMod(newMesh);
meshMod.changeMesh(mesh, false);
}
pointPriority = meshFilter.pointPriority();
}
if (collapseFaces)
{
meshFilterPtr.reset(new polyMeshFilter(mesh, pointPriority));
polyMeshFilter& meshFilter = meshFilterPtr();
const autoPtr<fvMesh>& newMesh = meshFilter.filteredMesh();
// Filter faces. Pass in the number of bad faces that are present
// from the previous edge filtering to use as a stopping criterion.
meshFilter.filter(nBadFaces);
@ -141,6 +185,8 @@ int main(int argc, char *argv[])
meshMod.changeMesh(mesh, false);
}
pointPriority = meshFilter.pointPriority();
}
// Write resulting mesh
@ -157,6 +203,7 @@ int main(int argc, char *argv[])
<< runTime.timeName() << nl << endl;
mesh.write();
pointPriority.write();
}
Info<< nl << "ExecutionTime = " << runTime.elapsedCpuTime() << " s"

9
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C
View File

@ -966,6 +966,15 @@ Foam::DistributedDelaunayMesh<Triangulation>::rangeInsertReferredWithInfo
) << "Point is outside affine hull! pt = " << pointToInsert
<< endl;
}
else if (lt == Triangulation::OUTSIDE_CONVEX_HULL)
{
// @todo Can this be optimised?
//
// Only want to insert if a connection is formed between
// pointToInsert and an internal or internal boundary point.
hint = Triangulation::insert(pointToInsert, c);
inserted = true;
}
else
{
// Get the cells that conflict with p in a vector V,

20
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
View File

@ -41,9 +41,26 @@ License
namespace Foam
{
defineTypeNameAndDebug(conformalVoronoiMesh, 0);
defineTypeNameAndDebug(conformalVoronoiMesh, 0);
template<>
const char* NamedEnum
<
conformalVoronoiMesh::dualMeshPointType,
5
>::names[] =
{
"internal",
"surface",
"featureEdge",
"featurePoint",
"constrained"
};
}
const Foam::NamedEnum<Foam::conformalVoronoiMesh::dualMeshPointType, 5>
Foam::conformalVoronoiMesh::dualMeshPointTypeNames_;
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
@ -1600,7 +1617,6 @@ void Foam::conformalVoronoiMesh::move()
printVertexInfo(Info);
}
// Write the intermediate mesh, do not filter the dual faces.
if (time().outputTime())
{
writeMesh(time().timeName());

17
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
View File

@ -115,6 +115,19 @@ public:
typedef List<pointIndexHitAndFeature> pointIndexHitAndFeatureList;
typedef DynamicList<pointIndexHitAndFeature> pointIndexHitAndFeatureDynList;
// Static data
enum dualMeshPointType
{
internal = 0,
surface = 1,
featureEdge = 2,
featurePoint = 3,
constrained = 4
};
static const NamedEnum<dualMeshPointType, 5> dualMeshPointTypeNames_;
private:
@ -682,14 +695,12 @@ private:
//- Merge vertices that are identical
void mergeIdenticalDualVertices
(
const pointField& pts,
const labelList& boundaryPts
const pointField& pts
);
label mergeIdenticalDualVertices
(
const pointField& pts,
const labelList& boundaryPts,
Map<label>& dualPtIndexMap
) const;

460
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
View File

@ -207,209 +207,216 @@ void Foam::conformalVoronoiMesh::checkCells()
}
void Foam::conformalVoronoiMesh::checkDuals()
{
List<List<Point> > pointFieldList(Pstream::nProcs());
List<Point> duals(number_of_finite_cells());
// PackedBoolList bPoints(number_of_finite_cells());
// indexDualVertices(duals, bPoints);
label count = 0;//duals.size();
duals.setSize(number_of_finite_cells());
globalIndex gIndex(number_of_vertices());
for
(
Delaunay::Finite_cells_iterator cit = finite_cells_begin();
cit != finite_cells_end();
++cit
)
{
if (cit->hasFarPoint())
{
continue;
}
duals[count++] = cit->circumcenter();
// List<labelPair> cellVerticesPair(4);
// List<Point> cellVertices(4);
//void Foam::conformalVoronoiMesh::checkDuals()
//{
// List<List<Point> > pointFieldList(Pstream::nProcs());
//
// for (label vI = 0; vI < 4; ++vI)
// List<Point> duals(number_of_finite_cells());
//
// typedef CGAL::Exact_predicates_exact_constructions_kernel EK2;
// typedef CGAL::Regular_triangulation_euclidean_traits_3<EK2> EK;
// typedef CGAL::Cartesian_converter<baseK::Kernel, EK2> To_exact;
// typedef CGAL::Cartesian_converter<EK2, baseK::Kernel> Back_from_exact;
//
//// PackedBoolList bPoints(number_of_finite_cells());
//
//// indexDualVertices(duals, bPoints);
//
// label count = 0;//duals.size();
//
// duals.setSize(number_of_finite_cells());
//
// globalIndex gIndex(number_of_vertices());
//
// for
// (
// Delaunay::Finite_cells_iterator cit = finite_cells_begin();
// cit != finite_cells_end();
// ++cit
// )
// {
// if (cit->hasFarPoint())
// {
// cellVerticesPair[vI] = labelPair
// (
// cit->vertex(vI)->procIndex(),
// cit->vertex(vI)->index()
// );
// cellVertices[vI] = cit->vertex(vI)->point();
// continue;
// }
//
// labelList oldToNew;
// sortedOrder(cellVerticesPair, oldToNew);
// oldToNew = invert(oldToNew.size(), oldToNew);
// inplaceReorder(oldToNew, cellVerticesPair);
// inplaceReorder(oldToNew, cellVertices);
// duals[count++] = cit->circumcenter();
//
// duals[count++] = CGAL::circumcenter
// (
// cellVertices[0],
// cellVertices[1],
// cellVertices[2],
// cellVertices[3]
// );
// To_exact to_exact;
// Back_from_exact back_from_exact;
// EK::Construct_circumcenter_3 exact_circumcenter =
// EK().construct_circumcenter_3_object();
//// List<labelPair> cellVerticesPair(4);
//// List<Point> cellVertices(4);
////
//// for (label vI = 0; vI < 4; ++vI)
//// {
//// cellVerticesPair[vI] = labelPair
//// (
//// cit->vertex(vI)->procIndex(),
//// cit->vertex(vI)->index()
//// );
//// cellVertices[vI] = cit->vertex(vI)->point();
//// }
////
//// labelList oldToNew;
//// sortedOrder(cellVerticesPair, oldToNew);
//// oldToNew = invert(oldToNew.size(), oldToNew);
//// inplaceReorder(oldToNew, cellVerticesPair);
//// inplaceReorder(oldToNew, cellVertices);
////
//// duals[count++] = CGAL::circumcenter
//// (
//// cellVertices[0],
//// cellVertices[1],
//// cellVertices[2],
//// cellVertices[3]
//// );
//
// duals[count++] = topoint
// (
// back_from_exact
// (
// exact_circumcenter
// (
// to_exact(cit->vertex(0)->point()),
// to_exact(cit->vertex(1)->point()),
// to_exact(cit->vertex(2)->point()),
// to_exact(cit->vertex(3)->point())
// )
// )
// );
}
Pout<< "Duals Calculated " << count << endl;
duals.setSize(count);
pointFieldList[Pstream::myProcNo()] = duals;
Pstream::gatherList(pointFieldList);
if (Pstream::master())
{
Info<< "Checking on master processor the dual locations of each " << nl
<< "processor point list against the master dual list." << nl
<< "There are " << pointFieldList.size() << " processors" << nl
<< "The size of each processor's dual list is:" << endl;
forAll(pointFieldList, pfI)
{
Info<< " Proc " << pfI << " has " << pointFieldList[pfI].size()
<< " duals" << endl;
}
label nNonMatches = 0;
label nNearMatches = 0;
label nExactMatches = 0;
forAll(pointFieldList[0], pI)
{
const Point& masterPoint = pointFieldList[0][pI];
bool foundMatch = false;
bool foundNearMatch = false;
scalar minCloseness = GREAT;
Point closestPoint(0, 0, 0);
forAll(pointFieldList, pfI)
{
if (pfI == 0)
{
continue;
}
// label pfI = 1;
forAll(pointFieldList[pfI], pISlave)
{
const Point& slavePoint
= pointFieldList[pfI][pISlave];
if (masterPoint == slavePoint)
{
foundMatch = true;
break;
}
const scalar closeness = mag
(
topoint(masterPoint) - topoint(slavePoint)
);
if (closeness < 1e-12)
{
foundNearMatch = true;
}
else
{
if (closeness < minCloseness)
{
minCloseness = closeness;
closestPoint = slavePoint;
}
}
}
if (!foundMatch)
{
if (foundNearMatch)
{
CGAL::Gmpq x(CGAL::to_double(masterPoint.x()));
CGAL::Gmpq y(CGAL::to_double(masterPoint.y()));
CGAL::Gmpq z(CGAL::to_double(masterPoint.z()));
std::cout<< "master = " << x << " " << y << " " << z
<< std::endl;
CGAL::Gmpq xs(CGAL::to_double(closestPoint.x()));
CGAL::Gmpq ys(CGAL::to_double(closestPoint.y()));
CGAL::Gmpq zs(CGAL::to_double(closestPoint.z()));
std::cout<< "slave = " << xs << " " << ys << " " << zs
<< std::endl;
nNearMatches++;
}
else
{
nNonMatches++;
Info<< " Closest point to " << masterPoint << " is "
<< closestPoint << nl
<< " Separation is " << minCloseness << endl;
CGAL::Gmpq x(CGAL::to_double(masterPoint.x()));
CGAL::Gmpq y(CGAL::to_double(masterPoint.y()));
CGAL::Gmpq z(CGAL::to_double(masterPoint.z()));
std::cout<< "master = " << x << " " << y << " " << z
<< std::endl;
CGAL::Gmpq xs(CGAL::to_double(closestPoint.x()));
CGAL::Gmpq ys(CGAL::to_double(closestPoint.y()));
CGAL::Gmpq zs(CGAL::to_double(closestPoint.z()));
std::cout<< "slave = " << xs << " " << ys << " " << zs
<< std::endl;
}
}
else
{
nExactMatches++;
}
}
}
Info<< "Found " << nNonMatches << " non-matching duals" << nl
<< " and " << nNearMatches << " near matches"
<< " and " << nExactMatches << " exact matches" << endl;
}
}
//// To_exact to_exact;
//// Back_from_exact back_from_exact;
//// EK::Construct_circumcenter_3 exact_circumcenter =
//// EK().construct_circumcenter_3_object();
////
//// duals[count++] = topoint
//// (
//// back_from_exact
//// (
//// exact_circumcenter
//// (
//// to_exact(cit->vertex(0)->point()),
//// to_exact(cit->vertex(1)->point()),
//// to_exact(cit->vertex(2)->point()),
//// to_exact(cit->vertex(3)->point())
//// )
//// )
//// );
// }
//
// Pout<< "Duals Calculated " << count << endl;
//
// duals.setSize(count);
//
// pointFieldList[Pstream::myProcNo()] = duals;
//
// Pstream::gatherList(pointFieldList);
//
// if (Pstream::master())
// {
// Info<< "Checking on master processor the dual locations of each" << nl
// << " processor point list against the master dual list." << nl
// << "There are " << pointFieldList.size() << " processors" << nl
// << "The size of each processor's dual list is:" << endl;
//
// forAll(pointFieldList, pfI)
// {
// Info<< " Proc " << pfI << " has " << pointFieldList[pfI].size()
// << " duals" << endl;
// }
//
// label nNonMatches = 0;
// label nNearMatches = 0;
// label nExactMatches = 0;
//
// forAll(pointFieldList[0], pI)
// {
// const Point& masterPoint = pointFieldList[0][pI];
//
// bool foundMatch = false;
// bool foundNearMatch = false;
//
// scalar minCloseness = GREAT;
// Point closestPoint(0, 0, 0);
//
// forAll(pointFieldList, pfI)
// {
// if (pfI == 0)
// {
// continue;
// }
//
//// label pfI = 1;
//
// forAll(pointFieldList[pfI], pISlave)
// {
// const Point& slavePoint
// = pointFieldList[pfI][pISlave];
//
// if (masterPoint == slavePoint)
// {
// foundMatch = true;
// break;
// }
//
// const scalar closeness = mag
// (
// topoint(masterPoint) - topoint(slavePoint)
// );
//
// if (closeness < 1e-12)
// {
// foundNearMatch = true;
// }
// else
// {
// if (closeness < minCloseness)
// {
// minCloseness = closeness;
// closestPoint = slavePoint;
// }
// }
// }
//
// if (!foundMatch)
// {
// if (foundNearMatch)
// {
// CGAL::Gmpq x(CGAL::to_double(masterPoint.x()));
// CGAL::Gmpq y(CGAL::to_double(masterPoint.y()));
// CGAL::Gmpq z(CGAL::to_double(masterPoint.z()));
//
// std::cout<< "master = " << x << " " << y << " " << z
// << std::endl;
//
// CGAL::Gmpq xs(CGAL::to_double(closestPoint.x()));
// CGAL::Gmpq ys(CGAL::to_double(closestPoint.y()));
// CGAL::Gmpq zs(CGAL::to_double(closestPoint.z()));
// std::cout<< "slave = " << xs << " " << ys << " "
// << zs
// << std::endl;
//
// nNearMatches++;
// }
// else
// {
// nNonMatches++;
// Info<< "Closest point to " << masterPoint << " is "
// << closestPoint << nl
// << " Separation is " << minCloseness << endl;
//
// CGAL::Gmpq x(CGAL::to_double(masterPoint.x()));
// CGAL::Gmpq y(CGAL::to_double(masterPoint.y()));
// CGAL::Gmpq z(CGAL::to_double(masterPoint.z()));
//
// std::cout<< "master = " << x << " " << y << " " << z
// << std::endl;
//
// CGAL::Gmpq xs(CGAL::to_double(closestPoint.x()));
// CGAL::Gmpq ys(CGAL::to_double(closestPoint.y()));
// CGAL::Gmpq zs(CGAL::to_double(closestPoint.z()));
// std::cout<< "slave = " << xs << " " << ys << " "
// << zs
// << std::endl;
// }
// }
// else
// {
// nExactMatches++;
// }
// }
// }
//
// Info<< "Found " << nNonMatches << " non-matching duals" << nl
// << " and " << nNearMatches << " near matches"
// << " and " << nExactMatches << " exact matches" << endl;
// }
//}
void Foam::conformalVoronoiMesh::checkVertices()
@ -578,7 +585,7 @@ void Foam::conformalVoronoiMesh::calcDualMesh
Info<< nl << "Merging identical points" << endl;
// There is no guarantee that a merge of close points is no-risk
mergeIdenticalDualVertices(points, boundaryPts);
mergeIdenticalDualVertices(points);
}
// Final dual face and owner neighbour construction
@ -813,8 +820,7 @@ void Foam::conformalVoronoiMesh::calcTetMesh
void Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
(
const pointField& pts,
const labelList& boundaryPts
const pointField& pts
)
{
// Assess close points to be merged
@ -829,7 +835,6 @@ void Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
nPtsMerged = mergeIdenticalDualVertices
(
pts,
boundaryPts,
dualPtIndexMap
);
@ -851,7 +856,6 @@ void Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
Foam::label Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
(
const pointField& pts,
const labelList& boundaryPts,
Map<label>& dualPtIndexMap
) const
{
@ -883,6 +887,19 @@ Foam::label Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
if (p1 == p2)
{
// if (c1->parallelDualVertex() || c2->parallelDualVertex())
// {
// if (c1->vertexLowestProc() < c2->vertexLowestProc())
// {
// dualPtIndexMap.insert(c1I, c1I);
// dualPtIndexMap.insert(c2I, c1I);
// }
// else
// {
// dualPtIndexMap.insert(c1I, c2I);
// dualPtIndexMap.insert(c2I, c2I);
// }
// }
if (c1I < c2I)
{
dualPtIndexMap.insert(c1I, c1I);
@ -1338,13 +1355,13 @@ void Foam::conformalVoronoiMesh::checkCellSizing()
timeCheck("Start of Cell Sizing");
labelList boundaryPts(number_of_finite_cells(), -1);
labelList boundaryPts(number_of_finite_cells(), internal);
pointField ptsField;
indexDualVertices(ptsField, boundaryPts);
// Merge close dual vertices.
mergeIdenticalDualVertices(ptsField, boundaryPts);
mergeIdenticalDualVertices(ptsField);
autoPtr<polyMesh> meshPtr = createPolyMeshFromPoints(ptsField);
const polyMesh& pMesh = meshPtr();
@ -1755,7 +1772,7 @@ void Foam::conformalVoronoiMesh::indexDualVertices
boundaryPts.setSize
(
number_of_finite_cells() + nConstrainedVertices,
-1
internal
);
if (foamyHexMeshControls().guardFeaturePoints())
@ -1774,7 +1791,7 @@ void Foam::conformalVoronoiMesh::indexDualVertices
topoint(vit->point());
boundaryPts[number_of_finite_cells() + nConstrainedVertices] =
1;
constrained;
nConstrainedVertices++;
}
@ -1974,15 +1991,40 @@ void Foam::conformalVoronoiMesh::indexDualVertices
if (cit->boundaryDualVertex())
{
if (cit->featureEdgeDualVertex())
if (cit->featurePointDualVertex())
{
boundaryPts[cit->cellIndex()] = 1;
boundaryPts[cit->cellIndex()] = featurePoint;
}
else if (cit->featureEdgeDualVertex())
{
boundaryPts[cit->cellIndex()] = featureEdge;
}
else
{
boundaryPts[cit->cellIndex()] = 0;
boundaryPts[cit->cellIndex()] = surface;
}
}
else if
(
cit->baffleBoundaryDualVertex()
)
{
boundaryPts[cit->cellIndex()] = surface;
}
else if
(
cit->vertex(0)->featureEdgePoint()
&& cit->vertex(1)->featureEdgePoint()
&& cit->vertex(2)->featureEdgePoint()
&& cit->vertex(3)->featureEdgePoint()
)
{
boundaryPts[cit->cellIndex()] = featureEdge;
}
else
{
boundaryPts[cit->cellIndex()] = internal;
}
}
else
{

202
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C
View File

@ -958,7 +958,7 @@ void Foam::conformalVoronoiMesh::writeMesh
{
Info<< nl << "Filtering edges on polyMesh" << nl << endl;
meshFilter.reset(new polyMeshFilter(mesh));
meshFilter.reset(new polyMeshFilter(mesh, boundaryPts));
// Filter small edges only. This reduces the number of faces so that
// the face filtering is sped up.
@ -974,9 +974,28 @@ void Foam::conformalVoronoiMesh::writeMesh
if (foamyHexMeshControls().filterFaces())
{
labelIOList boundaryPtsIO
(
IOobject
(
"pointPriority",
instance,
time(),
IOobject::NO_READ,
IOobject::NO_WRITE
),
labelList(mesh.nPoints(), labelMin)
);
forAll(mesh.points(), ptI)
{
boundaryPtsIO[ptI] = mesh.pointZones().whichZone(ptI);
}
Info<< nl << "Filtering faces on polyMesh" << nl << endl;
meshFilter.reset(new polyMeshFilter(mesh));
meshFilter.reset(new polyMeshFilter(mesh, boundaryPtsIO));
meshFilter().filter(nInitialBadFaces);
{
@ -1005,159 +1024,44 @@ void Foam::conformalVoronoiMesh::writeMesh
<< endl;
}
// volTensorField alignments
// (
// IOobject
// (
// "alignmentsField",
// runTime_.timeName(),
// runTime_,
// IOobject::NO_READ,
// IOobject::AUTO_WRITE
// ),
// mesh,
// tensor::zero
// );
//
// forAll(mesh.cellCentres(), pI)
// {
// Vertex_handle nearV =
// nearest_vertex
// (
// toPoint<Point>(mesh.cellCentres()[pI])
// );
// alignments[pI] = nearV->alignment();
// }
// alignments.write();
//
// {
// volVectorField alignmentx
// (
// IOobject
// (
// "alignmentsx",
// runTime_.timeName(),
// runTime_,
// IOobject::NO_READ,
// IOobject::AUTO_WRITE
// ),
// mesh,
// vector::zero
// );
// forAll(alignmentx, aI)
// {
// alignmentx[aI] = alignments[aI].x();
// }
// alignmentx.write();
// }
// {
// volVectorField alignmenty
// (
// IOobject
// (
// "alignmentsy",
// runTime_.timeName(),
// runTime_,
// IOobject::NO_READ,
// IOobject::AUTO_WRITE
// ),
// mesh,
// vector::zero
// );
// forAll(alignmenty, aI)
// {
// alignmenty[aI] = alignments[aI].y();
// }
// alignmenty.write();
// }
// {
// volVectorField alignmentz
// (
// IOobject
// (
// "alignmentsz",
// runTime_.timeName(),
// runTime_,
// IOobject::NO_READ,
// IOobject::AUTO_WRITE
// ),
// mesh,
// vector::zero
// );
// forAll(alignmentz, aI)
// {
// alignmentz[aI] = alignments[aI].z();
// }
// alignmentz.write();
// }
labelIOList boundaryIOPts
(
IOobject
(
"boundaryPoints",
instance,
runTime_,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
boundaryPts
);
// Dump list of boundary points
forAll(mesh.boundaryMesh(), patchI)
{
const polyPatch& pp = mesh.boundaryMesh()[patchI];
pointScalarField boundaryPtsScalarField
(
IOobject
(
"boundaryPoints_collapsed",
instance,
time(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
pointMesh::New(mesh),
scalar(labelMin)
);
if (!isA<coupledPolyPatch>(pp))
labelIOList boundaryPtsIO
(
IOobject
(
"pointPriority",
instance,
time(),
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
labelList(mesh.nPoints(), labelMin)
);
forAll(mesh.points(), ptI)
{
forAll(pp, fI)
{
const face& boundaryFace = pp[fI];
forAll(boundaryFace, pI)
{
const label boundaryPointI = boundaryFace[pI];
boundaryIOPts[boundaryPointI] = boundaryPts[boundaryPointI];
}
}
boundaryPtsScalarField[ptI] = mesh.pointZones().whichZone(ptI);
boundaryPtsIO[ptI] = mesh.pointZones().whichZone(ptI);
}
boundaryPtsScalarField.write();
boundaryPtsIO.write();
}
boundaryIOPts.write();
// forAllConstIter(labelHashSet, pointsInPatch, pI)
// {
// const Foam::point& ptMaster = mesh.points()[pI.key()];
//
// forAllConstIter(labelHashSet, pointsInPatch, ptI)
// {
// if (ptI.key() != pI.key())
// {
// const Foam::point& ptSlave = mesh.points()[ptI.key()];
//
// const scalar dist = mag(ptMaster - ptSlave);
// if (ptMaster == ptSlave)
// {
// Pout<< "Point(" << pI.key() << ") " << ptMaster
// << " == "
// << "(" << ptI.key() << ") " << ptSlave
// << endl;
// }
// else if (dist == 0)
// {
// Pout<< "Point(" << pI.key() << ") " << ptMaster
// << " ~= "
// << "(" << ptI.key() << ") " << ptSlave
// << endl;
// }
// }
// }
// }
// writeCellSizes(mesh);
// writeCellAlignments(mesh);

5
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
View File

@ -176,6 +176,9 @@ public:
//- Does the Delaunay cell have a far point
inline bool hasFarPoint() const;
//- Does the Delaunay cell have a referred point
inline bool hasReferredPoint() const;
//- Does the Delaunay cell have a feature point
inline bool hasFeaturePoint() const;
@ -216,6 +219,8 @@ public:
// least one Delaunay vertex outside and at least one inside
inline bool boundaryDualVertex() const;
inline bool baffleBoundaryDualVertex() const;
//- A dual vertex on a feature edge will result from this Delaunay cell
inline bool featureEdgeDualVertex() const;

81
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
View File

@ -189,6 +189,19 @@ inline bool CGAL::indexedCell<Gt, Cb>::hasFarPoint() const
}
template<class Gt, class Cb>
inline bool CGAL::indexedCell<Gt, Cb>::hasReferredPoint() const
{
return
(
this->vertex(0)->referred()
|| this->vertex(1)->referred()
|| this->vertex(2)->referred()
|| this->vertex(3)->referred()
);
}
template<class Gt, class Cb>
inline bool CGAL::indexedCell<Gt, Cb>::hasFeaturePoint() const
{
@ -372,6 +385,14 @@ inline bool CGAL::indexedCell<Gt, Cb>::anyInternalOrBoundaryDualVertex() const
template<class Gt, class Cb>
inline bool CGAL::indexedCell<Gt, Cb>::boundaryDualVertex() const
{
// return
// (
// this->vertex(0)->boundaryPoint()
// && this->vertex(1)->boundaryPoint()
// && this->vertex(2)->boundaryPoint()
// && this->vertex(3)->boundaryPoint()
// );
return
(
(
@ -387,6 +408,41 @@ inline bool CGAL::indexedCell<Gt, Cb>::boundaryDualVertex() const
|| this->vertex(3)->externalBoundaryPoint()
)
);
// Foam::label nBoundaryPoints = 0;
//
// for (Foam::label i = 0; i < 4; ++i)
// {
// Vertex_handle v = this->vertex(i);
//
// if (v->boundaryPoint())
// {
// nBoundaryPoints++;
// }
// }
//
// return (nBoundaryPoints > 1);
}
template<class Gt, class Cb>
inline bool CGAL::indexedCell<Gt, Cb>::baffleBoundaryDualVertex() const
{
return
(
(
this->vertex(0)->internalBafflePoint()
|| this->vertex(1)->internalBafflePoint()
|| this->vertex(2)->internalBafflePoint()
|| this->vertex(3)->internalBafflePoint()
)
&& (
this->vertex(0)->externalBafflePoint()
|| this->vertex(1)->externalBafflePoint()
|| this->vertex(2)->externalBafflePoint()
|| this->vertex(3)->externalBafflePoint()
)
);
}
@ -400,6 +456,31 @@ inline bool CGAL::indexedCell<Gt, Cb>::featureEdgeDualVertex() const
&& this->vertex(2)->featureEdgePoint()
&& this->vertex(3)->featureEdgePoint()
);
// (
// this->vertex(0)->featureEdgePoint()
// || this->vertex(1)->featureEdgePoint()
// || this->vertex(2)->featureEdgePoint()
// || this->vertex(3)->featureEdgePoint()
// )
// && (
// (
// this->vertex(0)->featureEdgePoint()
// || this->vertex(0)->featurePoint()
// )
// && (
// this->vertex(1)->featureEdgePoint()
// || this->vertex(1)->featurePoint()
// )
// && (
// this->vertex(2)->featureEdgePoint()
// || this->vertex(2)->featurePoint()
// )
// && (
// this->vertex(3)->featureEdgePoint()
// || this->vertex(3)->featurePoint()
// )
// )
// );
}

2
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.H
View File

@ -243,8 +243,10 @@ public:
inline bool surfacePoint() const;
inline bool internalBoundaryPoint() const;
inline bool internalBafflePoint() const;
inline bool externalBoundaryPoint() const;
inline bool externalBafflePoint() const;
inline bool constrained() const;

8
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexEnum.C
View File

@ -30,13 +30,17 @@ License
template<>
const char*
Foam::NamedEnum<Foam::indexedVertexEnum::vertexType, 11>::names[] =
Foam::NamedEnum<Foam::indexedVertexEnum::vertexType, 15>::names[] =
{
"Unassigned",
"Internal",
"InternalNearBoundary",
"InternalSurface",
"InternalSurfaceBaffle",
"ExternalSurfaceBaffle",
"InternalFeatureEdge",
"InternalFeatureEdgeBaffle",
"ExternalFeatureEdgeBaffle",
"InternalFeaturePoint",
"ExternalSurface",
"ExternalFeatureEdge",
@ -45,7 +49,7 @@ Foam::NamedEnum<Foam::indexedVertexEnum::vertexType, 11>::names[] =
"Constrained"
};
const Foam::NamedEnum<Foam::indexedVertexEnum::vertexType, 11>
const Foam::NamedEnum<Foam::indexedVertexEnum::vertexType, 15>
Foam::indexedVertexEnum::vertexTypeNames_;

28
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexEnum.H
View File

@ -49,17 +49,21 @@ public:
enum vertexType
{
vtUnassigned = 0,
vtInternal = 1,
vtInternalNearBoundary = 2,
vtInternalSurface = 3,
vtInternalFeatureEdge = 4,
vtInternalFeaturePoint = 5,
vtExternalSurface = 6,
vtExternalFeatureEdge = 7,
vtExternalFeaturePoint = 8,
vtFar = 9,
vtConstrained = 10
vtUnassigned = 0,
vtInternal = 1,
vtInternalNearBoundary = 2,
vtInternalSurface = 3,
vtInternalSurfaceBaffle = 4,
vtExternalSurfaceBaffle = 5,
vtInternalFeatureEdge = 6,
vtInternalFeatureEdgeBaffle = 7,
vtExternalFeatureEdgeBaffle = 8,
vtInternalFeaturePoint = 9,
vtExternalSurface = 10,
vtExternalFeatureEdge = 11,
vtExternalFeaturePoint = 12,
vtFar = 13,
vtConstrained = 14
};
enum vertexMotion
@ -68,7 +72,7 @@ public:
movable = 1
};
static const Foam::NamedEnum<vertexType, 11> vertexTypeNames_;
static const Foam::NamedEnum<vertexType, 15> vertexTypeNames_;
static const Foam::NamedEnum<vertexMotion, 2> vertexMotionNames_;

20
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexI.H
View File

@ -307,6 +307,16 @@ inline bool CGAL::indexedVertex<Gt, Vb>::internalBoundaryPoint() const
return type_ >= vtInternalSurface && type_ <= vtInternalFeaturePoint;
}
template<class Gt, class Vb>
inline bool CGAL::indexedVertex<Gt, Vb>::internalBafflePoint() const
{
return
(
type_ == vtInternalSurfaceBaffle
|| type_ == vtInternalFeatureEdgeBaffle
);
}
template<class Gt, class Vb>
inline bool CGAL::indexedVertex<Gt, Vb>::externalBoundaryPoint() const
@ -314,6 +324,16 @@ inline bool CGAL::indexedVertex<Gt, Vb>::externalBoundaryPoint() const
return type_ >= vtExternalSurface && type_ <= vtExternalFeaturePoint;
}
template<class Gt, class Vb>
inline bool CGAL::indexedVertex<Gt, Vb>::externalBafflePoint() const
{
return
(
type_ == vtExternalSurfaceBaffle
|| type_ == vtExternalFeatureEdgeBaffle
);
}
template<class Gt, class Vb>
inline bool CGAL::indexedVertex<Gt, Vb>::constrained() const

3
applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/rayShooting/rayShooting.C
View File

@ -49,7 +49,6 @@ void rayShooting::splitLine
{
Foam::point midPoint(l.centre());
const scalar localCellSize(cellShapeControls().cellSize(midPoint));
const scalar lineLength(l.mag());
const scalar minDistFromSurfaceSqr
(
@ -64,6 +63,8 @@ void rayShooting::splitLine
)
{
// Add extra points if line length is much bigger than local cell size
// const scalar lineLength(l.mag());
//
// if (lineLength > 4.0*localCellSize)
// {
// splitLine

43
src/OpenFOAM/db/error/IOerror.C
View File

@ -255,28 +255,31 @@ void Foam::IOerror::abort()
Foam::Ostream& Foam::operator<<(Ostream& os, const IOerror& ioErr)
{
os << endl
<< ioErr.title().c_str() << endl
<< ioErr.message().c_str() << endl << endl;
os << "file: " << ioErr.ioFileName().c_str();
if (ioErr.ioStartLineNumber() >= 0 && ioErr.ioEndLineNumber() >= 0)
if (!os.bad())
{
os << " from line " << ioErr.ioStartLineNumber()
<< " to line " << ioErr.ioEndLineNumber() << '.';
}
else if (ioErr.ioStartLineNumber() >= 0)
{
os << " at line " << ioErr.ioStartLineNumber() << '.';
}
os << endl
<< ioErr.title().c_str() << endl
<< ioErr.message().c_str() << endl << endl;
if (IOerror::level >= 2 && ioErr.sourceFileLineNumber())
{
os << endl << endl
<< " From function " << ioErr.functionName().c_str() << endl
<< " in file " << ioErr.sourceFileName().c_str()
<< " at line " << ioErr.sourceFileLineNumber() << '.';
os << "file: " << ioErr.ioFileName().c_str();
if (ioErr.ioStartLineNumber() >= 0 && ioErr.ioEndLineNumber() >= 0)
{
os << " from line " << ioErr.ioStartLineNumber()
<< " to line " << ioErr.ioEndLineNumber() << '.';
}
else if (ioErr.ioStartLineNumber() >= 0)
{
os << " at line " << ioErr.ioStartLineNumber() << '.';
}
if (IOerror::level >= 2 && ioErr.sourceFileLineNumber())
{
os << endl << endl
<< " From function " << ioErr.functionName().c_str() << endl
<< " in file " << ioErr.sourceFileName().c_str()
<< " at line " << ioErr.sourceFileLineNumber() << '.';
}
}
return os;

6
src/OpenFOAM/meshes/treeBoundBox/treeBoundBox.C
View File

@ -269,7 +269,8 @@ bool Foam::treeBoundBox::intersects
const direction endBits = posBits(end);
pt = start;
while (true)
// Allow maximum of 3 clips.
for (label i = 0; i < 4; ++i)
{
direction ptBits = posBits(pt);
@ -380,6 +381,9 @@ bool Foam::treeBoundBox::intersects
}
}
}
// Can end up here if the end point is on the edge of the boundBox
return true;
}

1
src/dynamicMesh/Make/files
View File

@ -99,6 +99,7 @@ createShellMesh/createShellMesh.C
extrudePatchMesh/extrudePatchMesh.C
polyMeshFilter/polyMeshFilterSettings.C
polyMeshFilter/polyMeshFilter.C
LIB = $(FOAM_LIBBIN)/libdynamicMesh

1575
src/dynamicMesh/polyMeshFilter/polyMeshFilter.C
View File

File diff suppressed because it is too large Load Diff

85
src/dynamicMesh/polyMeshFilter/polyMeshFilter.H
View File

@ -43,6 +43,7 @@ SourceFiles
#include "List.H"
#include "autoPtr.H"
#include "scalarField.H"
#include "polyMeshFilterSettings.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -59,6 +60,8 @@ class faceZone;
\*---------------------------------------------------------------------------*/
class polyMeshFilter
:
private polyMeshFilterSettings
{
// Private data
@ -68,52 +71,9 @@ class polyMeshFilter
//- Copy of the original mesh to perform the filtering on
autoPtr<fvMesh> newMeshPtr_;
//- Dictionary containing the coefficient sub-dictionaries
const IOdictionary dict_;
//- After collapsing, check the mesh quality and redo the collapsing
// iteration if there are too many bad faces in the mesh
Switch controlMeshQuality_;
//- Coefficients for collapsing edges
const dictionary& collapseEdgesCoeffDict_;
//- Coefficients for collapsing faces
dictionary collapseFacesCoeffDict_;
//- Coefficients for controlling the mesh quality
dictionary meshQualityCoeffDict_;
//- Remove edges shorter than this length
const scalar minLen_;
//- Merge points that are only attached to two edges and have an angle
// between the edge greater than this value
const scalar maxCos_;
//- The amount that the local minimum edge length will be reduced by if
// the edge is part of a collapse string that generates poor quality
// faces
const scalar edgeReductionFactor_;
//- Maximum number of outer iterations
const label maxIterations_;
//- Maximum number of smoothing iterations of minEdgeLen_ and
// faceFilterFactor_
const label maxSmoothIters_;
//- Initialisation value of faceFilterFactor_
const scalar initialFaceLengthFactor_;
//- The amount that the local face size factor will be reduced by if
// the face is part of a collapse string that generates poor quality
// faces
const scalar faceReductionFactor_;
//- Maximum number of times a deleted point can be associated with the
// creation of a bad face it is forced to be kept.
const label maxPointErrorCount_;
//-
labelList originalPointPriority_;
autoPtr<labelList> pointPriority_;
//- The minimum edge length for each edge
scalarField minEdgeLen_;
@ -124,6 +84,22 @@ class polyMeshFilter
// Private Member Functions
label filterFacesLoop(const label nOriginalBadFaces);
label filterFaces
(
polyMesh& newMesh,
scalarField& newMeshFaceFilterFactor,
labelList& origToCurrentPointMap
);
label filterEdges
(
polyMesh& newMesh,
scalarField& newMeshMinEdgeLen,
labelList& origToCurrentPointMap
);
//- Increment pointErrorCount for points attached to a bad face
void updatePointErrorCount
(
@ -160,7 +136,11 @@ class polyMeshFilter
// + >0 : point on a processor patch with that ID
// @todo Need to mark boundaryEdges as well, as an edge may have two
// boundary points but not itself lie on a boundary
labelList findBoundaryPoints(const polyMesh& mesh) const;
void updatePointPriorities
(
const polyMesh& newMesh,
const labelList& pointMap
);
//- Print min/mean/max data for a field
void printScalarFieldStats
@ -213,6 +193,9 @@ public:
//- Construct from fvMesh
explicit polyMeshFilter(const fvMesh& mesh);
//- Construct from fvMesh and a label list of point priorities
polyMeshFilter(const fvMesh& mesh, const labelList& pointPriority);
//- Destructor
~polyMeshFilter();
@ -226,6 +209,8 @@ public:
// mesh has actually been filtered.
const autoPtr<fvMesh>& filteredMesh() const;
const autoPtr<labelList>& pointPriority() const;
// Edit
@ -235,11 +220,11 @@ public:
//- Filter edges and faces
label filter(const label nOriginalBadFaces);
//- Filter all faces that are in the face zone indirectPatchFaces
label filter(const faceZone& fZone);
//- Filter edges only.
label filterEdges(const label nOriginalBadFaces);
//- Filter all faces that are in the face zone indirectPatchFaces
label filterFaceZone(const faceZone& fZone);
};

276
src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C
View File

@ -182,53 +182,46 @@ void Foam::edgeCollapser::collapseToEdge
Map<point>& collapsePointToLocation
) const
{
const face& f = mesh_.faces()[faceI];
// Negative half
Foam::point collapseToPtA =
collapseAxis*(sum(dNeg)/dNeg.size() - dShift) + fC;
Foam::point collapseToPtA(GREAT, GREAT, GREAT);
//collapseAxis*(sum(dNeg)/dNeg.size() - dShift) + fC;
DynamicList<label> faceBoundaryPts(f.size());
DynamicList<label> faceFeaturePts(f.size());
label maxPriority = labelMin;
DynamicList<label> maxPriorityPts(max(dNeg.size(), dPos.size()));
forAll(facePtsNeg, fPtI)
{
if (pointPriority[facePtsNeg[fPtI]] == 1)
const label facePointI = facePtsNeg[fPtI];
const label facePtPriority = pointPriority[facePointI];
if (facePtPriority > maxPriority)
{
faceFeaturePts.append(facePtsNeg[fPtI]);
maxPriority = facePtPriority;
maxPriorityPts.clear();
maxPriorityPts.append(facePointI);
}
else if (pointPriority[facePtsNeg[fPtI]] == 0)
else if (facePtPriority == maxPriority)
{
faceBoundaryPts.append(facePtsNeg[fPtI]);
maxPriorityPts.append(facePointI);
}
}
if (!faceBoundaryPts.empty() || !faceFeaturePts.empty())
if (!maxPriorityPts.empty())
{
if (!faceFeaturePts.empty())
{
collapseToPtA = pts[faceFeaturePts.first()];
}
else if (faceBoundaryPts.size() == 2)
{
collapseToPtA =
0.5
*(
pts[faceBoundaryPts[0]]
+ pts[faceBoundaryPts[1]]
);
}
else if (faceBoundaryPts.size() <= f.size())
{
face bFace(faceBoundaryPts);
Foam::point averagePt(vector::zero);
collapseToPtA = bFace.centre(pts);
forAll(maxPriorityPts, ptI)
{
averagePt += pts[maxPriorityPts[ptI]];
}
collapseToPtA = averagePt/maxPriorityPts.size();
// collapseToPtA = pts[maxPriorityPts.first()];
}
faceFeaturePts.clear();
faceBoundaryPts.clear();
maxPriority = labelMin;
maxPriorityPts.clear();
labelList faceEdgesNeg = edgesFromPoints(faceI, facePtsNeg);
@ -244,47 +237,37 @@ void Foam::edgeCollapser::collapseToEdge
// Positive half
Foam::point collapseToPtB
= collapseAxis*(sum(dPos)/dPos.size() - dShift) + fC;
Foam::point collapseToPtB(GREAT, GREAT, GREAT);
// = collapseAxis*(sum(dPos)/dPos.size() - dShift) + fC;
forAll(facePtsPos, fPtI)
{
if (pointPriority[facePtsPos[fPtI]] == 1)
const label facePointI = facePtsPos[fPtI];
const label facePtPriority = pointPriority[facePointI];
if (facePtPriority > maxPriority)
{
faceFeaturePts.append(facePtsPos[fPtI]);
maxPriority = facePtPriority;
maxPriorityPts.clear();
maxPriorityPts.append(facePointI);
}
else if (pointPriority[facePtsPos[fPtI]] == 0)
else if (facePtPriority == maxPriority)
{
// If there is a point which is on the boundary,
// use it as the point to collapse others to, will
// use the first boundary point encountered if
// there are multiple boundary points.
faceBoundaryPts.append(facePtsPos[fPtI]);
maxPriorityPts.append(facePointI);
}
}
if (!faceBoundaryPts.empty() || !faceFeaturePts.empty())
if (!maxPriorityPts.empty())
{
if (!faceFeaturePts.empty())
{
collapseToPtB = pts[faceFeaturePts.first()];
}
else if (faceBoundaryPts.size() == 2)
{
collapseToPtB =
0.5
*(
pts[faceBoundaryPts[0]]
+ pts[faceBoundaryPts[1]]
);
}
else if (faceBoundaryPts.size() <= f.size())
{
face bFace(faceBoundaryPts);
Foam::point averagePt(vector::zero);
collapseToPtB = bFace.centre(pts);
forAll(maxPriorityPts, ptI)
{
averagePt += pts[maxPriorityPts[ptI]];
}
collapseToPtB = averagePt/maxPriorityPts.size();
// collapseToPtB = pts[maxPriorityPts.first()];
}
labelList faceEdgesPos = edgesFromPoints(faceI, facePtsPos);
@ -316,45 +299,79 @@ void Foam::edgeCollapser::collapseToPoint
Foam::point collapseToPt = fC;
DynamicList<label> faceBoundaryPts(f.size());
DynamicList<label> faceFeaturePts(f.size());
label maxPriority = labelMin;
DynamicList<label> maxPriorityPts(f.size());
forAll(facePts, fPtI)
{
if (pointPriority[facePts[fPtI]] == 1)
const label facePointI = facePts[fPtI];
const label facePtPriority = pointPriority[facePointI];
if (facePtPriority > maxPriority)
{
faceFeaturePts.append(facePts[fPtI]);
maxPriority = facePtPriority;
maxPriorityPts.clear();
maxPriorityPts.append(facePointI);
}
else if (pointPriority[facePts[fPtI]] == 0)
else if (facePtPriority == maxPriority)
{
faceBoundaryPts.append(facePts[fPtI]);
maxPriorityPts.append(facePointI);
}
}
if (!faceBoundaryPts.empty() || !faceFeaturePts.empty())
if (!maxPriorityPts.empty())
{
if (!faceFeaturePts.empty())
{
collapseToPt = pts[faceFeaturePts.first()];
}
else if (faceBoundaryPts.size() == 2)
{
collapseToPt =
0.5
*(
pts[faceBoundaryPts[0]]
+ pts[faceBoundaryPts[1]]
);
}
else if (faceBoundaryPts.size() <= f.size())
{
face bFace(faceBoundaryPts);
Foam::point averagePt(vector::zero);
collapseToPt = bFace.centre(pts);
forAll(maxPriorityPts, ptI)
{
averagePt += pts[maxPriorityPts[ptI]];
}
collapseToPt = averagePt/maxPriorityPts.size();
// collapseToPt = pts[maxPriorityPts.first()];
}
const labelList faceEdges = mesh_.faceEdges()[faceI];
// DynamicList<label> faceBoundaryPts(f.size());
// DynamicList<label> faceFeaturePts(f.size());
//
// forAll(facePts, fPtI)
// {
// if (pointPriority[facePts[fPtI]] == 1)
// {
// faceFeaturePts.append(facePts[fPtI]);
// }
// else if (pointPriority[facePts[fPtI]] == 0)
// {
// faceBoundaryPts.append(facePts[fPtI]);
// }
// }
//
// if (!faceBoundaryPts.empty() || !faceFeaturePts.empty())
// {
// if (!faceFeaturePts.empty())
// {
// collapseToPt = pts[faceFeaturePts.first()];
// }
// else if (faceBoundaryPts.size() == 2)
// {
// collapseToPt =
// 0.5
// *(
// pts[faceBoundaryPts[0]]
// + pts[faceBoundaryPts[1]]
// );
// }
// else if (faceBoundaryPts.size() <= f.size())
// {
// face bFace(faceBoundaryPts);
//
// collapseToPt = bFace.centre(pts);
// }
// }
const labelList& faceEdges = mesh_.faceEdges()[faceI];
forAll(faceEdges, eI)
{
@ -740,35 +757,51 @@ Foam::label Foam::edgeCollapser::edgeMaster
{
label masterPoint = -1;
label e0 = e.start();
label e1 = e.end();
const label e0 = e.start();
const label e1 = e.end();
// Collapse edge to point with higher priority.
if (pointPriority[e0] >= 0)
const label e0Priority = pointPriority[e0];
const label e1Priority = pointPriority[e1];
if (e0Priority > e1Priority)
{
if (pointPriority[e1] >= 0)
{
// Both points have high priority. Choose one to collapse to.
// Note: should look at feature edges/points!
masterPoint = e0;
}
else
{
masterPoint = e0;
}
masterPoint = e0;
}
else
else if (e0Priority < e1Priority)
{
if (pointPriority[e1] >= 0)
{
masterPoint = e1;
}
else
{
// None on boundary. Neither is a master.
return -1;
}
masterPoint = e1;
}
else if (e0Priority == e1Priority)
{
masterPoint = e0;
}
// // Collapse edge to point with higher priority.
// if (pointPriority[e0] >= 0)
// {
// if (pointPriority[e1] >= 0)
// {
// // Both points have high priority. Choose one to collapse to.
// // Note: should look at feature edges/points!
// masterPoint = e0;
// }
// else
// {
// masterPoint = e0;
// }
// }
// else
// {
// if (pointPriority[e1] >= 0)
// {
// masterPoint = e1;
// }
// else
// {
// // None on boundary. Neither is a master.
// return -1;
// }
// }
return masterPoint;
}
@ -784,7 +817,10 @@ void Foam::edgeCollapser::checkBoundaryPointMergeEdges
{
const pointField& points = mesh_.points();
if (pointPriority[pointI] >= 0 && pointPriority[otherPointI] < 0)
const label e0Priority = pointPriority[pointI];
const label e1Priority = pointPriority[otherPointI];
if (e0Priority > e1Priority)
{
collapsePointToLocation.set
(
@ -792,7 +828,7 @@ void Foam::edgeCollapser::checkBoundaryPointMergeEdges
points[pointI]
);
}
else
else if (e0Priority < e1Priority)
{
collapsePointToLocation.set
(
@ -800,6 +836,22 @@ void Foam::edgeCollapser::checkBoundaryPointMergeEdges
points[otherPointI]
);
}
else // e0Priority == e1Priority
{
collapsePointToLocation.set
(
pointI,
points[otherPointI]
);
// Foam::point averagePt
// (
// 0.5*(points[otherPointI] + points[pointI])
// );
//
// collapsePointToLocation.set(pointI, averagePt);
// collapsePointToLocation.set(otherPointI, averagePt);
}
}
@ -1963,7 +2015,7 @@ Foam::labelPair Foam::edgeCollapser::markSmallSliverFaces
{
const face& f = faces[fI];
if (faceFilterFactor[fI] == 0)
if (faceFilterFactor[fI] <= 0)
{
continue;
}
@ -2030,7 +2082,7 @@ Foam::labelPair Foam::edgeCollapser::markFaceZoneEdges
const face& f = faces[fI];
if (faceFilterFactor[fI] == 0)
if (faceFilterFactor[fI] <= 0)
{
continue;
}