diff --git a/applications/test/PtrList/Test-PtrList.C b/applications/test/PtrList/Test-PtrList.C
index 9241460a59d12c2c23439fae1178d07ad0378fd8..664dd5b501ef35dca4b3bbd063eb080b9eaf3431 100644
--- a/applications/test/PtrList/Test-PtrList.C
+++ b/applications/test/PtrList/Test-PtrList.C
@@ -32,6 +32,8 @@ Description
#include "scalar.H"
#include "IOstreams.H"
#include "PtrList.H"
+#include "plane.H"
+#include "DynamicList.H"
using namespace Foam;
@@ -76,6 +78,7 @@ int main(int argc, char *argv[])
{
PtrList<Scalar> list1(10);
PtrList<Scalar> list2(15);
+ PtrList<Scalar> listApp;
forAll(list1, i)
{
@@ -87,9 +90,14 @@ int main(int argc, char *argv[])
list2.set(i, new Scalar(10 + 1.3*i));
}
+ for (label i = 0; i < 5; ++i)
+ {
+ listApp.append(new Scalar(1.3*i));
+ }
Info<<"list1: " << list1 << endl;
Info<<"list2: " << list2 << endl;
+ Info<<"listApp: " << listApp << endl;
Info<<"indirectly delete some items via set(.., 0) :" << endl;
for (label i = 0; i < 3; i++)
@@ -115,6 +123,13 @@ int main(int argc, char *argv[])
Info<<"list2: " << list2 << endl;
Info<<"list3: " << list3 << endl;
+ PtrList<plane> planes;
+ planes.append(new plane(vector::one, vector::one));
+ planes.append(new plane(vector(1,2,3), vector::one));
+
+ forAll(planes, p)
+ Info<< "plane " << planes[p] << endl;
+
Info<< nl << "Done." << endl;
return 0;
}
diff --git a/applications/test/sort/Test-sortList.C b/applications/test/sort/Test-sortList.C
index 0521159ebbf3f2d672fbecbc5d64dc24286e55a5..d92a4a7ad429ba9d1ce16fb0759676e6ac8d24d5 100644
--- a/applications/test/sort/Test-sortList.C
+++ b/applications/test/sort/Test-sortList.C
@@ -50,11 +50,18 @@ int main(int argc, char *argv[])
labelList a(orig);
sortedOrder(a, order);
+ SortableList<label> aReverse(a.size());
+ aReverse = a;
+
Info<< "unsorted: " << a << endl;
sort(a);
Info<< "sorted: " << a << endl;
Info<< "indices: " << order << endl;
+ aReverse.reverseSort();
+ Info<< "reverse sorted: " << aReverse << endl;
+ Info<< "reverse indices: " << aReverse.indices() << endl;
+
SortableList<label> b(orig);
Info<< "unsorted: " << orig << endl;
Info<< "sorted: " << b << endl;
diff --git a/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C b/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
index 0515205c7887c4cc27fe36f8c359335049fa63d2..5b24606613cc87af2a4a1d24ffe378163ea84ebb 100644
--- a/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
+++ b/applications/utilities/mesh/advanced/collapseEdges/collapseEdges.C
@@ -52,6 +52,7 @@ Usage
#include "polyTopoChange.H"
#include "fvMesh.H"
#include "polyMeshFilter.H"
+#include "faceSet.H"
using namespace Foam;
@@ -74,9 +75,9 @@ int main(int argc, char *argv[])
argList::addOption
(
- "collapseFaceZone",
- "zoneName",
- "Collapse faces that are in the supplied face zone"
+ "collapseFaceSet",
+ "faceSet",
+ "Collapse faces that are in the supplied face set"
);
# include "addOverwriteOption.H"
@@ -93,15 +94,15 @@ int main(int argc, char *argv[])
const bool overwrite = args.optionFound("overwrite");
const bool collapseFaces = args.optionFound("collapseFaces");
- const bool collapseFaceZone = args.optionFound("collapseFaceZone");
+ const bool collapseFaceSet = args.optionFound("collapseFaceSet");
- if (collapseFaces && collapseFaceZone)
+ if (collapseFaces && collapseFaceSet)
{
FatalErrorIn("main(int, char*[])")
<< "Both face zone collapsing and face collapsing have been"
<< "selected. Choose only one of:" << nl
<< " -collapseFaces" << nl
- << " -collapseFaceZone <faceZoneName>"
+ << " -collapseFaceSet <faceSet>"
<< abort(FatalError);
}
@@ -117,7 +118,6 @@ int main(int argc, char *argv[])
),
labelList(mesh.nPoints(), labelMin)
);
-
forAll(timeDirs, timeI)
{
runTime.setTime(timeDirs[timeI], timeI);
@@ -128,6 +128,14 @@ int main(int argc, char *argv[])
label nBadFaces = 0;
+ faceSet indirectPatchFaces
+ (
+ mesh,
+ "indirectPatchFaces",
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ );
+
{
meshFilterPtr.set(new polyMeshFilter(mesh, pointPriority));
polyMeshFilter& meshFilter = meshFilterPtr();
@@ -139,19 +147,26 @@ int main(int argc, char *argv[])
// the face filtering is sped up.
nBadFaces = meshFilter.filterEdges(0);
{
- polyTopoChange meshMod(newMesh);
+ polyTopoChange meshMod(newMesh());
meshMod.changeMesh(mesh, false);
+
+ polyMeshFilter::copySets(newMesh(), mesh);
}
pointPriority = meshFilter.pointPriority();
}
- if (collapseFaceZone)
+ if (collapseFaceSet)
{
- const word faceZoneName = args.optionRead<word>("collapseFaceZone");
-
- const faceZone& fZone = mesh.faceZones()[faceZoneName];
+ const word faceSetName(args.optionRead<word>("collapseFaceSet"));
+ faceSet fSet
+ (
+ mesh,
+ faceSetName,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ );
meshFilterPtr.reset(new polyMeshFilter(mesh, pointPriority));
polyMeshFilter& meshFilter = meshFilterPtr();
@@ -160,11 +175,13 @@ int main(int argc, char *argv[])
// 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(fZone);
+ meshFilter.filter(fSet);
{
polyTopoChange meshMod(newMesh);
meshMod.changeMesh(mesh, false);
+
+ polyMeshFilter::copySets(newMesh(), mesh);
}
pointPriority = meshFilter.pointPriority();
@@ -184,6 +201,8 @@ int main(int argc, char *argv[])
polyTopoChange meshMod(newMesh);
meshMod.changeMesh(mesh, false);
+
+ polyMeshFilter::copySets(newMesh(), mesh);
}
pointPriority = meshFilter.pointPriority();
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/Make/options b/applications/utilities/mesh/generation/foamyHexMesh/Make/options
index f03dbf90976d5f7ecedf098b24beced6e5cb968c..256fa40164d063935dd9cbfc23553679b0d59b46 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyHexMesh/Make/options
@@ -10,8 +10,10 @@ include $(GENERAL_RULES)/CGAL
EXE_INC = \
${EXE_FROUNDING_MATH} \
${EXE_NDEBUG} \
+ ${CGAL_EXACT} \
${CGAL_INEXACT} \
${CGAL_INC} \
+ ${c++CGALWARN} \
-IconformalVoronoiMesh/lnInclude \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
index 3afdf9cf139312594e7fbd6ad83069c3d59824b4..8f90bb2118e4f04709e69d280dfb6530f57fb25c 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.C
@@ -190,20 +190,25 @@ void Foam::DelaunayMesh<Triangulation>::reset()
resetVertexCount();
resetCellCount();
- insertPoints(vertices);
+ insertPoints(vertices, false);
Info<< "Inserted " << vertexCount() << " fixed points" << endl;
}
template<class Triangulation>
-void Foam::DelaunayMesh<Triangulation>::insertPoints(const List<Vb>& vertices)
+Foam::Map<Foam::label> Foam::DelaunayMesh<Triangulation>::insertPoints
+(
+ const List<Vb>& vertices,
+ const bool reIndex
+)
{
- rangeInsertWithInfo
+ return rangeInsertWithInfo
(
vertices.begin(),
vertices.end(),
- false
+ false,
+ reIndex
);
}
@@ -268,7 +273,7 @@ const
template<class Triangulation>
template<class PointIterator>
-void Foam::DelaunayMesh<Triangulation>::rangeInsertWithInfo
+Foam::Map<Foam::label> Foam::DelaunayMesh<Triangulation>::rangeInsertWithInfo
(
PointIterator begin,
PointIterator end,
@@ -307,6 +312,10 @@ void Foam::DelaunayMesh<Triangulation>::rangeInsertWithInfo
Vertex_handle hint;
+ Map<label> oldToNewIndex(points.size());
+
+ label maxIndex = -1;
+
for
(
typename vectorPairPointIndex::const_iterator p = points.begin();
@@ -335,11 +344,14 @@ void Foam::DelaunayMesh<Triangulation>::rangeInsertWithInfo
{
if (reIndex)
{
+ const label oldIndex = vert.index();
hint->index() = getNewVertexIndex();
+ oldToNewIndex.insert(oldIndex, hint->index());
}
else
{
hint->index() = vert.index();
+ maxIndex = max(maxIndex, vert.index());
}
hint->type() = vert.type();
hint->procIndex() = vert.procIndex();
@@ -347,6 +359,13 @@ void Foam::DelaunayMesh<Triangulation>::rangeInsertWithInfo
hint->alignment() = vert.alignment();
}
}
+
+ if (!reIndex)
+ {
+ vertexCount_ = maxIndex + 1;
+ }
+
+ return oldToNewIndex;
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
index a7cbf26766e48b08436610c67802973a284eb1d2..0bf843f1c9074bd61a457db9bac1d4991b4024f2 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DelaunayMesh.H
@@ -228,7 +228,11 @@ public:
void reset();
//- Insert the list of vertices (calls rangeInsertWithInfo)
- void insertPoints(const List<Vb>& vertices);
+ Map<label> insertPoints
+ (
+ const List<Vb>& vertices,
+ const bool reIndex
+ );
//- Function inserting points into a triangulation and setting the
// index and type data of the point in the correct order. This is
@@ -237,7 +241,7 @@ public:
// Adapted from a post on the CGAL lists: 2010-01/msg00004.html by
// Sebastien Loriot (Geometry Factory).
template<class PointIterator>
- void rangeInsertWithInfo
+ Map<label> rangeInsertWithInfo
(
PointIterator begin,
PointIterator end,
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C
index 40c76aec9daf18cb5ba766a3aaa36490f827127e..b598a7f5de7ae72435dd0045680faa7a9c086a03 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/DelaunayMesh/DistributedDelaunayMesh.C
@@ -536,7 +536,8 @@ Foam::label Foam::DistributedDelaunayMesh<Triangulation>::referVertices
labelPairHashSet pointsNotInserted = rangeInsertReferredWithInfo
(
referredVertices.begin(),
- referredVertices.end()
+ referredVertices.end(),
+ true
);
if (!pointsNotInserted.empty())
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/files b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/files
index efeef3586643c824c81a7363a3aae6f056948a92..190b77775b07a563d071120f1e4853515cfe863a 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/files
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/files
@@ -76,6 +76,7 @@ faceAreaWeightModel/piecewiseLinearRamp/piecewiseLinearRamp.C
searchableSurfaceFeatures/searchableSurfaceFeatures.C
searchableSurfaceFeatures/searchableBoxFeatures.C
+searchableSurfaceFeatures/searchablePlateFeatures.C
searchableSurfaceFeatures/triSurfaceMeshFeatures.C
LIB = $(FOAM_LIBBIN)/libconformalVoronoiMesh
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/options b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/options
index 2be5855323b10652c1de0a434a8d1ec95f6f2c16..ac0dfdc4f8944e187dbc558adc076b5c44b0e203 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/Make/options
@@ -11,8 +11,10 @@ FFLAGS = -DCGAL_FILES='"${CGAL_ARCH_PATH}/share/files"'
EXE_INC = \
${EXE_FROUNDING_MATH} \
${EXE_NDEBUG} \
+ ${CGAL_EXACT} \
${CGAL_INEXACT} \
${CGAL_INC} \
+ ${c++CGALWARN} \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/parallel/decompose/decompositionMethods/lnInclude \
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
index 71ba6933fbb3a67780bafad1113238cd8b4e3b00..2bfe84404ad84a35171757ae95b74e3f868f0257 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.C
@@ -628,49 +628,6 @@ Foam::tensorField Foam::cellShapeControlMesh::dumpAlignments() const
}
-void Foam::cellShapeControlMesh::insertBoundingPoints
-(
- const boundBox& bb,
- const cellSizeAndAlignmentControls& sizeControls
-)
-{
- // Loop over bound box points and get cell size and alignment
- const pointField bbPoints(bb.points());
-
- forAll(bbPoints, pI)
- {
- const Foam::point& pt = bbPoints[pI];
-
- // Cell size here will return default cell size
- const scalar cellSize = sizeControls.cellSize(pt);
-
- if (debug)
- {
- Info<< "Insert Bounding Point: " << pt << " " << cellSize << endl;
- }
-
- // Get the cell size of the nearest surface.
-// geometryToConformTo_.findSurfaceNearest
-// (
-// pt,
-// GREAT,
-// surfHit,
-// hitSurface
-// );
-
- const tensor alignment = tensor::I;
-
- insert
- (
- pt,
- cellSize,
- alignment,
- Vb::vtInternalNearBoundary
- );
- }
-}
-
-
void Foam::cellShapeControlMesh::write() const
{
Info<< "Writing " << meshSubDir << endl;
@@ -781,4 +738,78 @@ void Foam::cellShapeControlMesh::write() const
}
+Foam::label Foam::cellShapeControlMesh::estimateCellCount
+(
+ const autoPtr<backgroundMeshDecomposition>& decomposition
+) const
+{
+ // Loop over all the tets and estimate the cell count in each one
+
+ scalar cellCount = 0;
+
+ for
+ (
+ Finite_cells_iterator cit = finite_cells_begin();
+ cit != finite_cells_end();
+ ++cit
+ )
+ {
+ if (!cit->hasFarPoint() && !is_infinite(cit))
+ {
+ // @todo Check if tet centre is on the processor..
+ CGAL::Tetrahedron_3<baseK> tet
+ (
+ cit->vertex(0)->point(),
+ cit->vertex(1)->point(),
+ cit->vertex(2)->point(),
+ cit->vertex(3)->point()
+ );
+
+ pointFromPoint centre = topoint(CGAL::centroid(tet));
+
+ if
+ (
+ Pstream::parRun()
+ && !decomposition().positionOnThisProcessor(centre)
+ )
+ {
+ continue;
+ }
+
+ scalar volume = CGAL::to_double(tet.volume());
+
+ scalar averagedPointCellSize = 0;
+ //scalar averagedPointCellSize = 1;
+
+ // Get an average volume by averaging the cell size of the vertices
+ for (label vI = 0; vI < 4; ++vI)
+ {
+ averagedPointCellSize += cit->vertex(vI)->targetCellSize();
+ //averagedPointCellSize *= cit->vertex(vI)->targetCellSize();
+ }
+
+ averagedPointCellSize /= 4;
+ //averagedPointCellSize = ::sqrt(averagedPointCellSize);
+
+// if (averagedPointCellSize < SMALL)
+// {
+// Pout<< "Volume = " << volume << endl;
+//
+// for (label vI = 0; vI < 4; ++vI)
+// {
+// Pout<< "Point " << vI
+// << ", point = " << topoint(cit->vertex(vI)->point())
+// << ", size = " << cit->vertex(vI)->targetCellSize()
+// << endl;
+// }
+// }
+
+ cellCount += volume/pow(averagedPointCellSize, 3);
+ }
+ }
+
+ return cellCount;
+}
+
+
// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
index e6af5bde0f6520065df2c8b264cb1aa3461e71a1..79a397c49a3f75beafb3a6d5d6f9512c1c964d6e 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellShapeControlMesh/cellShapeControlMesh.H
@@ -153,15 +153,14 @@ public:
tensorField dumpAlignments() const;
- void insertBoundingPoints
- (
- const boundBox& bb,
- const cellSizeAndAlignmentControls& sizeControls
- );
-
void writeTriangulation();
void write() const;
+
+ label estimateCellCount
+ (
+ const autoPtr<backgroundMeshDecomposition>& decomposition
+ ) const;
};
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
index bff35cff69a19ff1202ba881e63e645a236db52e..49fc665e934c7c4686a34996b06c7c7ce47c0e9c 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/cellSizeAndAlignmentControl/searchableSurfaceControl/searchableSurfaceControl.C
@@ -429,6 +429,11 @@ void Foam::searchableSurfaceControl::initialVertices
vectorField normals(1);
searchableSurface_.getNormal(infoList, normals);
+ if (mag(normals[0]) < SMALL)
+ {
+ normals[0] = vector(1, 1, 1);
+ }
+
pointAlignment.set(new triad(normals[0]));
if (infoList[0].hit())
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
index 8d954799a1e26c49d21f96ba12091e8bf85b99bb..2be7d916d3589988d0b161ec2046767a2b1c7d3e 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/cellShapeControl/controlMeshRefinement/controlMeshRefinement.C
@@ -768,7 +768,7 @@ Foam::label Foam::controlMeshRefinement::refineMesh
}
}
- mesh_.insertPoints(verts);
+ mesh_.insertPoints(verts, false);
return verts.size();
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
index 10ea65f77213a0603ee7fd4f151a7ae6a4b1dd7b..c74740d81267723ece9c7e85cab6d4ef164cec71 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C
@@ -179,15 +179,23 @@ void Foam::conformalVoronoiMesh::insertInternalPoints
}
-void Foam::conformalVoronoiMesh::insertPoints
+Foam::Map<Foam::label> Foam::conformalVoronoiMesh::insertPointPairs
(
List<Vb>& vertices,
- bool distribute
+ bool distribute,
+ bool reIndex
)
{
if (Pstream::parRun() && distribute)
{
- decomposition_().distributePoints(vertices);
+ autoPtr<mapDistribute> mapDist =
+ decomposition_().distributePoints(vertices);
+
+ // Re-index the point pairs if one or both have been distributed.
+ // If both, remove
+
+ // If added a point, then need to know its point pair
+ // If one moved, then update procIndex locally
forAll(vertices, vI)
{
@@ -197,7 +205,8 @@ void Foam::conformalVoronoiMesh::insertPoints
label preReinsertionSize(number_of_vertices());
- this->DelaunayMesh<Delaunay>::insertPoints(vertices);
+ Map<label> oldToNewIndices =
+ this->DelaunayMesh<Delaunay>::insertPoints(vertices, reIndex);
const label nReinserted = returnReduce
(
@@ -208,17 +217,18 @@ void Foam::conformalVoronoiMesh::insertPoints
Info<< " Reinserted " << nReinserted << " vertices out of "
<< returnReduce(vertices.size(), sumOp<label>())
<< endl;
+
+ return oldToNewIndices;
}
void Foam::conformalVoronoiMesh::insertSurfacePointPairs
(
const pointIndexHitAndFeatureList& surfaceHits,
- const fileName fName
+ const fileName fName,
+ DynamicList<Vb>& pts
)
{
- DynamicList<Vb> pts(2.0*surfaceHits.size());
-
forAll(surfaceHits, i)
{
vectorField norm(1);
@@ -246,6 +256,7 @@ void Foam::conformalVoronoiMesh::insertSurfacePointPairs
pointPairDistance(surfacePt),
surfacePt,
normal,
+ true,
pts
);
}
@@ -256,6 +267,7 @@ void Foam::conformalVoronoiMesh::insertSurfacePointPairs
pointPairDistance(surfacePt),
surfacePt,
normal,
+ true,
pts
);
}
@@ -266,6 +278,7 @@ void Foam::conformalVoronoiMesh::insertSurfacePointPairs
pointPairDistance(surfacePt),
surfacePt,
-normal,
+ true,
pts
);
}
@@ -280,8 +293,6 @@ void Foam::conformalVoronoiMesh::insertSurfacePointPairs
}
}
- insertPoints(pts, true);
-
if (foamyHexMeshControls().objOutput() && fName != fileName::null)
{
DelaunayMeshTools::writeOBJ(time().path()/fName, pts);
@@ -292,11 +303,10 @@ void Foam::conformalVoronoiMesh::insertSurfacePointPairs
void Foam::conformalVoronoiMesh::insertEdgePointGroups
(
const pointIndexHitAndFeatureList& edgeHits,
- const fileName fName
+ const fileName fName,
+ DynamicList<Vb>& pts
)
{
- DynamicList<Vb> pts(3.0*edgeHits.size());
-
forAll(edgeHits, i)
{
if (edgeHits[i].first().hit())
@@ -318,10 +328,6 @@ void Foam::conformalVoronoiMesh::insertEdgePointGroups
}
}
- pts.shrink();
-
- insertPoints(pts, true);
-
if (foamyHexMeshControls().objOutput() && fName != fileName::null)
{
DelaunayMeshTools::writeOBJ(time().path()/fName, pts);
@@ -348,6 +354,28 @@ bool Foam::conformalVoronoiMesh::nearFeaturePt(const Foam::point& pt) const
}
+bool Foam::conformalVoronoiMesh::surfacePtNearFeatureEdge
+(
+ const Foam::point& pt
+) const
+{
+ scalar exclusionRangeSqr = surfacePtExclusionDistanceSqr(pt);
+
+ pointIndexHit info;
+ label featureHit;
+
+ geometryToConformTo_.findEdgeNearest
+ (
+ pt,
+ exclusionRangeSqr,
+ info,
+ featureHit
+ );
+
+ return info.hit();
+}
+
+
void Foam::conformalVoronoiMesh::insertInitialPoints()
{
Info<< nl << "Inserting initial points" << endl;
@@ -859,6 +887,7 @@ Foam::conformalVoronoiMesh::conformalVoronoiMesh
geometryToConformTo_
),
limitBounds_(),
+ ptPairs_(*this),
ftPtConformer_(*this),
edgeLocationTreePtr_(),
surfacePtLocationTreePtr_(),
@@ -1108,7 +1137,7 @@ void Foam::conformalVoronoiMesh::move()
vB->alignment().T() & cartesianDirections
);
- Field<vector> alignmentDirs(3);
+ Field<vector> alignmentDirs(alignmentDirsA);
forAll(alignmentDirsA, aA)
{
@@ -1266,7 +1295,7 @@ void Foam::conformalVoronoiMesh::move()
if
(
(
- (vA->internalPoint() || vB->internalPoint())
+ (vA->internalPoint() && vB->internalPoint())
&& (!vA->referred() || !vB->referred())
// ||
// (
@@ -1333,12 +1362,17 @@ void Foam::conformalVoronoiMesh::move()
// Point removal
if
(
- vA->internalPoint()
- && !vA->referred()
- && !vA->fixed()
- && vB->internalPoint()
- && !vB->referred()
- && !vB->fixed()
+ (
+ vA->internalPoint()
+ && !vA->referred()
+ && !vA->fixed()
+ )
+ &&
+ (
+ vB->internalPoint()
+ && !vB->referred()
+ && !vB->fixed()
+ )
)
{
// Only insert a point at the midpoint of
@@ -1507,9 +1541,7 @@ void Foam::conformalVoronoiMesh::move()
Info<< "Writing point displacement vectors to file." << endl;
OFstream str
(
- time().path()
- + "displacements_" + runTime_.timeName()
- + ".obj"
+ time().path()/"displacements_" + runTime_.timeName() + ".obj"
);
for
@@ -1607,6 +1639,74 @@ void Foam::conformalVoronoiMesh::move()
time().path()/"boundaryPoints_" + time().timeName() + ".obj",
*this
);
+
+ DelaunayMeshTools::writeOBJ
+ (
+ time().path()/"internalBoundaryPoints_" + time().timeName()
+ + ".obj",
+ *this,
+ Foam::indexedVertexEnum::vtInternalSurface,
+ Foam::indexedVertexEnum::vtInternalFeaturePoint
+ );
+
+ DelaunayMeshTools::writeOBJ
+ (
+ time().path()/"externalBoundaryPoints_" + time().timeName()
+ + ".obj",
+ *this,
+ Foam::indexedVertexEnum::vtExternalSurface,
+ Foam::indexedVertexEnum::vtExternalFeaturePoint
+ );
+
+ OBJstream multipleIntersections
+ (
+ "multipleIntersections_"
+ + time().timeName()
+ + ".obj"
+ );
+
+ for
+ (
+ Delaunay::Finite_edges_iterator eit = finite_edges_begin();
+ eit != finite_edges_end();
+ ++eit
+ )
+ {
+ Cell_handle c = eit->first;
+ Vertex_handle vA = c->vertex(eit->second);
+ Vertex_handle vB = c->vertex(eit->third);
+
+ Foam::point ptA(topoint(vA->point()));
+ Foam::point ptB(topoint(vB->point()));
+
+ List<pointIndexHit> surfHits;
+ labelList hitSurfaces;
+
+ geometryToConformTo().findSurfaceAllIntersections
+ (
+ ptA,
+ ptB,
+ surfHits,
+ hitSurfaces
+ );
+
+ if
+ (
+ surfHits.size() >= 2
+ || (
+ surfHits.size() == 0
+ && (
+ (vA->externalBoundaryPoint()
+ && vB->internalBoundaryPoint())
+ || (vB->externalBoundaryPoint()
+ && vA->internalBoundaryPoint())
+ )
+ )
+ )
+ {
+ multipleIntersections.write(linePointRef(ptA, ptB));
+ }
+ }
}
}
@@ -1629,49 +1729,6 @@ void Foam::conformalVoronoiMesh::move()
}
-//Foam::labelListList Foam::conformalVoronoiMesh::overlapsProc
-//(
-// const List<Foam::point>& centres,
-// const List<scalar>& radiusSqrs
-//) const
-//{
-// if (!Pstream::parRun())
-// {
-// return labelListList(centres.size(), labelList(0));
-// }
-//
-//// DynamicList<Foam::point> pts(number_of_vertices());
-//
-//// for
-//// (
-//// Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-//// vit != finite_vertices_end();
-//// vit++
-//// )
-//// {
-//// pts.append(topoint(vit->point()));
-//// }
-////
-//// dynamicIndexedOctree<dynamicTreeDataPoint> vertexOctree
-//// (
-//// dynamicTreeDataPoint(pts),
-//// treeBoundBox(min(pts), max(pts)),
-//// 10, // maxLevel
-//// 10, // leafSize
-//// 3.0 // duplicity
-//// );
-//
-// return decomposition_().overlapsProcessors
-// (
-// centres,
-// radiusSqrs,
-// *this,
-// false//,
-//// vertexOctree
-// );
-//}
-
-
void Foam::conformalVoronoiMesh::checkCoPlanarCells() const
{
typedef CGAL::Exact_predicates_exact_constructions_kernel Kexact;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
index 334330f970ed3fdb50d379325601a50bcb264275..2996493dbec4fa6f31fc9ca08a3c3d1f11ddb033 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H
@@ -77,7 +77,7 @@ SourceFiles
#include "Pair.H"
#include "DistributedDelaunayMesh.H"
#include "featurePointConformer.H"
-
+#include "pointPairs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -86,12 +86,10 @@ namespace Foam
// Forward declaration of classes
class initialPointsMethod;
-
class relaxationModel;
-
class faceAreaWeightModel;
-
class backgroundMeshDecomposition;
+class OBJstream;
/*---------------------------------------------------------------------------*\
Class conformalVoronoiMesh Declaration
@@ -165,6 +163,9 @@ private:
//- Limiting bound box before infinity begins
treeBoundBox limitBounds_;
+ //-
+ mutable pointPairs<Delaunay> ptPairs_;
+
//-
featurePointConformer ftPtConformer_;
@@ -229,10 +230,11 @@ private:
const bool distribute = false
);
- void insertPoints
+ Map<label> insertPointPairs
(
List<Vb>& vertices,
- bool distribute
+ bool distribute,
+ bool reIndex
);
//- Create a point-pair at a ppDist distance either side of
@@ -242,6 +244,7 @@ private:
const scalar ppDist,
const Foam::point& surfPt,
const vector& n,
+ const bool ptPair,
DynamicList<Vb>& pts
) const;
@@ -254,24 +257,20 @@ private:
const scalar ppDist,
const Foam::point& surfPt,
const vector& n,
+ const bool ptPair,
DynamicList<Vb>& pts
) const;
//- Check internal point is completely inside the meshable region
inline bool internalPointIsInside(const Foam::point& pt) const;
- inline bool isPointPair
- (
- const Vertex_handle& vA,
- const Vertex_handle& vB
- ) const;
-
//- Insert pairs of points on the surface with the given normals, at the
// specified spacing
void insertSurfacePointPairs
(
const pointIndexHitAndFeatureList& surfaceHits,
- const fileName fName = fileName::null
+ const fileName fName,
+ DynamicList<Vb>& pts
);
//- Insert groups of points to conform to an edge given a list of
@@ -280,7 +279,8 @@ private:
void insertEdgePointGroups
(
const pointIndexHitAndFeatureList& edgeHits,
- const fileName fName = fileName::null
+ const fileName fName,
+ DynamicList<Vb>& pts
);
void createEdgePointGroupByCirculating
@@ -351,6 +351,9 @@ private:
//- Check if a location is in exclusion range around a feature point
bool nearFeaturePt(const Foam::point& pt) const;
+ //- Check if a surface point is in exclusion range around a feature edge
+ bool surfacePtNearFeatureEdge(const Foam::point& pt) const;
+
//- Insert the initial points into the triangulation, based on the
// initialPointsMethod
void insertInitialPoints();
@@ -474,14 +477,11 @@ private:
label& hitSurface
) const;
- //- Find the "worst" incursion of the dual cell of a non-internal or
- // boundary point through the surface, subject to the
- // maxSurfaceProtrusion tolerance
void dualCellLargestSurfaceIncursion
(
const Delaunay::Finite_vertices_iterator& vit,
- pointIndexHit& surfHitLargest,
- label& hitSurfaceLargest
+ pointIndexHit& surfHit,
+ label& hitSurface
) const;
//- Write out vertex-processor occupancy information for debugging
@@ -695,7 +695,8 @@ private:
//- Merge vertices that are identical
void mergeIdenticalDualVertices
(
- const pointField& pts
+ const pointField& pts,
+ labelList& boundaryPts
);
label mergeIdenticalDualVertices
@@ -727,6 +728,8 @@ private:
// elements damage the mesh quality, allowing backtracking.
labelHashSet checkPolyMeshQuality(const pointField& pts) const;
+ label classifyBoundaryPoint(Cell_handle cit) const;
+
//- Index all of the the Delaunay cells and calculate their
//- dual points
void indexDualVertices
@@ -736,7 +739,11 @@ private:
);
//- Re-index all of the Delaunay cells
- void reindexDualVertices(const Map<label>& dualPtIndexMap);
+ void reindexDualVertices
+ (
+ const Map<label>& dualPtIndexMap,
+ labelList& boundaryPts
+ );
label createPatchInfo
(
@@ -846,6 +853,8 @@ private:
const PtrList<dictionary>& patchDicts
) const;
+ void writePointPairs(const fileName& fName) const;
+
//- Disallow default bitwise copy construct
conformalVoronoiMesh(const conformalVoronoiMesh&);
@@ -992,7 +1001,7 @@ public:
const wordList& patchNames,
const PtrList<dictionary>& patchDicts,
const pointField& cellCentres,
- const PackedBoolList& boundaryFacesToRemove
+ PackedBoolList& boundaryFacesToRemove
) const;
//- Calculate and write a field of the target cell size,
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
index cb17a1367dbb9410906f9a939c35a1ab3c734761..df2c8c2f651602d08c030867ce9916319b4bf69b 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C
@@ -30,487 +30,11 @@ License
#include "indexedCellChecks.H"
#include "OBJstream.H"
#include "indexedCellOps.H"
+#include "ListOps.H"
#include "DelaunayMeshTools.H"
-#include "CGAL/Exact_predicates_exact_constructions_kernel.h"
-#include "CGAL/Gmpq.h"
-
// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
-void Foam::conformalVoronoiMesh::checkCells()
-{
- List<List<FixedList<Foam::point, 4> > > cellListList(Pstream::nProcs());
-
- List<FixedList<Foam::point, 4> > cells(number_of_finite_cells());
-
- globalIndex gIndex(number_of_vertices());
-
- label count = 0;
- for
- (
- Delaunay::Finite_cells_iterator cit = finite_cells_begin();
- cit != finite_cells_end();
- ++cit
- )
- {
- if (tetrahedron(cit).volume() == 0)
- {
- Pout<< "ZERO VOLUME TET" << endl;
- Pout<< cit->info();
- Pout<< cit->dual();
- }
-
- if (cit->hasFarPoint())
- {
- continue;
- }
-
- List<labelPair> cellVerticesPair(4);
- List<Foam::point> cellVertices(4);
-
- for (label vI = 0; vI < 4; ++vI)
- {
- cellVerticesPair[vI] = labelPair
- (
- cit->vertex(vI)->procIndex(),
- cit->vertex(vI)->index()
- );
- cellVertices[vI] = topoint(cit->vertex(vI)->point());
- }
-
- List<Foam::point> cellVerticesOld(cellVertices);
- labelList oldToNew;
- sortedOrder(cellVerticesPair, oldToNew);
- oldToNew = invert(oldToNew.size(), oldToNew);
- inplaceReorder(oldToNew, cellVerticesPair);
- inplaceReorder(oldToNew, cellVertices);
-
-// FixedList<label, 4> globalTetCell
-// (
-// cit->globallyOrderedCellVertices(gIndex)
-// );
-//
-// FixedList<Point, 4> cellVertices(Point(0,0,0));
-//
-// forAll(globalTetCell, gvI)
-// {
-// label gI = globalTetCell[gvI];
-//
-// cellVertices[gvI] = cit->vertex(gI)->point();
-// }
-
-// if (cit->hasFarPoint())
-// {
-// continue;
-// }
-
- for (label i = 0; i < 4; ++i)
- {
- //cells[count][i] = topoint(cit->vertex(i)->point());
- cells[count][i] = cellVertices[i];
- }
-
- count++;
- }
-
- cells.setSize(count);
-
- cellListList[Pstream::myProcNo()] = cells;
-
- Pstream::gatherList(cellListList);
-
- if (Pstream::master())
- {
- Info<< "Checking on master processor the cells of each " << nl
- << "processor point list against the master cell list." << nl
- << "There are " << cellListList.size() << " processors" << nl
- << "The size of each processor's cell list is:" << endl;
-
- forAll(cellListList, cfI)
- {
- Info<< " Proc " << cfI << " has " << cellListList[cfI].size()
- << " cells" << endl;
- }
-
- label nMatches = 0, nMatchFoundDiffOrder = 0;
-
- forAll(cellListList[0], cmI)
- {
- const FixedList<Foam::point, 4>& masterCell = cellListList[0][cmI];
-
- bool matchFound = false;
- bool matchFoundDiffOrder = false;
-
- forAll(cellListList, cpI)
- {
- if (cpI == 0)
- {
- continue;
- }
-
- forAll(cellListList[cpI], csI)
- {
- const FixedList<Foam::point, 4>& slaveCell
- = cellListList[cpI][csI];
-
- if (masterCell == slaveCell)
- {
- matchFound = true;
- break;
- }
- else
- {
- label samePt = 0;
-
- forAll(masterCell, mI)
- {
- const Foam::point& mPt = masterCell[mI];
-
- forAll(slaveCell, sI)
- {
- const Foam::point& sPt = slaveCell[sI];
-
- if (mPt == sPt)
- {
- samePt++;
- }
- }
- }
-
- if (samePt == 4)
- {
- matchFoundDiffOrder = true;
-
- Pout<< masterCell << nl << slaveCell << endl;
-
- break;
- }
- }
- }
- }
-
- if (matchFound)
- {
- nMatches++;
- }
-
- if (matchFoundDiffOrder)
- {
- nMatchFoundDiffOrder++;
- }
- }
-
- Info<< "Found " << nMatches << " matching cells and "
- << nMatchFoundDiffOrder << " matching cells with different "
- << "vertex ordering"<< endl;
- }
-}
-
-
-//void Foam::conformalVoronoiMesh::checkDuals()
-//{
-// List<List<Point> > pointFieldList(Pstream::nProcs());
-//
-// 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())
-// {
-// continue;
-// }
-//
-// duals[count++] = cit->circumcenter();
-//
-//// 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]
-//// );
-//
-//// 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()
-{
- List<pointField> pointFieldList(Pstream::nProcs());
-
- pointField points(number_of_vertices());
-
- labelPairHashSet duplicateVertices;
-
- label count = 0;
- for
- (
- Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
- vit != finite_vertices_end();
- ++vit
- )
- {
- if (duplicateVertices.found(labelPair(vit->procIndex(), vit->index())))
- {
- Pout<< "DUPLICATE " << vit->procIndex() << vit->index() << endl;
- }
- else
- {
- duplicateVertices.insert(labelPair(vit->procIndex(), vit->index()));
- }
-
- points[count++] = topoint(vit->point());
- }
-
- pointFieldList[Pstream::myProcNo()] = points;
-
- Pstream::gatherList(pointFieldList);
-
- OFstream str("missingPoints.obj");
-
- if (Pstream::master())
- {
- Info<< "Checking on master processor the point locations of each " << nl
- << "processor point list against the master point list." << nl
- << "There are " << pointFieldList.size() << " processors" << nl
- << "The size of each processor's point list is:" << endl;
-
- forAll(pointFieldList, pfI)
- {
- Info<< " Proc " << pfI << " has " << pointFieldList[pfI].size()
- << " points" << endl;
- }
-
- label nNonMatches = 0;
-
- forAll(pointFieldList[0], pI)
- {
- const Foam::point& masterPoint = pointFieldList[0][pI];
-
- forAll(pointFieldList, pfI)
- {
- if (pI == 0)
- {
- continue;
- }
-
- bool foundMatch = false;
-
- forAll(pointFieldList[pfI], pISlave)
- {
- const Foam::point& slavePoint
- = pointFieldList[pfI][pISlave];
-
- if (masterPoint == slavePoint)
- {
- foundMatch = true;
- break;
- }
- }
-
- if (!foundMatch)
- {
- Info<< " Proc " << pfI << " Master != Slave -> "
- << masterPoint << endl;
-
- meshTools::writeOBJ(str, masterPoint);
-
- nNonMatches++;
- }
- }
- }
-
- Info<< "Found a total of " << nNonMatches << " non-matching points"
- << endl;
- }
-}
-
-
void Foam::conformalVoronoiMesh::calcDualMesh
(
pointField& points,
@@ -528,53 +52,6 @@ void Foam::conformalVoronoiMesh::calcDualMesh
{
timeCheck("Start calcDualMesh");
-// if (debug)
-// {
-// Pout<< nl << "Perfoming some checks . . ." << nl << nl
-// << "Total number of vertices = " << number_of_vertices() << nl
-// << "Total number of cells = " << number_of_finite_cells()
-// << endl;
-//
-// checkVertices();
-// checkCells();
-// checkDuals();
-//
-// Info<< nl << "Finished checks" << nl << endl;
-// }
-
-// OFstream str("attachedToFeature.obj");
-// label offset = 0;
-//
-// for
-// (
-// Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
-// vit != finite_vertices_end();
-// ++vit
-// )
-// {
-// if (vit->featurePoint())
-// {
-// std::list<Cell_handle> adjacentCells;
-//
-// finite_incident_cells(vit, std::back_inserter(adjacentCells));
-//
-// for
-// (
-// std::list<Cell_handle>::iterator acit = adjacentCells.begin();
-// acit != adjacentCells.end();
-// ++acit
-// )
-// {
-// if ((*acit)->real())
-// {
-// drawDelaunayCell(str, (*acit), offset);
-// offset++;
-//// meshTools::writeOBJ(str, topoint((*acit)->dual()));
-// }
-// }
-// }
-// }
-
setVertexSizeAndAlignment();
timeCheck("After setVertexSizeAndAlignment");
@@ -585,7 +62,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);
+ mergeIdenticalDualVertices(points, boundaryPts);
}
// Final dual face and owner neighbour construction
@@ -820,7 +297,8 @@ void Foam::conformalVoronoiMesh::calcTetMesh
void Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
(
- const pointField& pts
+ const pointField& pts,
+ labelList& boundaryPts
)
{
// Assess close points to be merged
@@ -838,7 +316,7 @@ void Foam::conformalVoronoiMesh::mergeIdenticalDualVertices
dualPtIndexMap
);
- reindexDualVertices(dualPtIndexMap);
+ reindexDualVertices(dualPtIndexMap, boundaryPts);
reduce(nPtsMerged, sumOp<label>());
@@ -1245,7 +723,6 @@ Foam::conformalVoronoiMesh::createPolyMeshFromPoints
false
);
- //createCellCentres(cellCentres);
cellCentres = DelaunayMeshTools::allPoints(*this);
labelList cellToDelaunayVertex(removeUnusedCells(owner, neighbour));
@@ -1327,24 +804,6 @@ Foam::conformalVoronoiMesh::createPolyMeshFromPoints
pMesh.addPatches(patches);
- // Info<< "ADDPATCHES NOT IN PARALLEL" << endl;
-
- // forAll(patches, p)
- // {
- // patches[p] = new polyPatch
- // (
- // patchNames[p],
- // patchSizes[p],
- // patchStarts[p],
- // p,
- // pMesh.boundaryMesh()
- // );
- // }
-
- // pMesh.addPatches(patches, false);
-
- // pMesh.overrideCellCentres(cellCentres);
-
return meshPtr;
}
@@ -1361,7 +820,7 @@ void Foam::conformalVoronoiMesh::checkCellSizing()
indexDualVertices(ptsField, boundaryPts);
// Merge close dual vertices.
- mergeIdenticalDualVertices(ptsField);
+ mergeIdenticalDualVertices(ptsField, boundaryPts);
autoPtr<polyMesh> meshPtr = createPolyMeshFromPoints(ptsField);
const polyMesh& pMesh = meshPtr();
@@ -1740,6 +1199,41 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::checkPolyMeshQuality
}
+Foam::label Foam::conformalVoronoiMesh::classifyBoundaryPoint
+(
+ Cell_handle cit
+) const
+{
+ if (cit->boundaryDualVertex())
+ {
+ if (cit->featurePointDualVertex())
+ {
+ return featurePoint;
+ }
+ else if (cit->featureEdgeDualVertex())
+ {
+ return featureEdge;
+ }
+ else
+ {
+ return surface;
+ }
+ }
+ else if (cit->baffleSurfaceDualVertex())
+ {
+ return surface;
+ }
+ else if (cit->baffleEdgeDualVertex())
+ {
+ return featureEdge;
+ }
+ else
+ {
+ return internal;
+ }
+}
+
+
void Foam::conformalVoronoiMesh::indexDualVertices
(
pointField& pts,
@@ -1989,42 +1483,7 @@ void Foam::conformalVoronoiMesh::indexDualVertices
// }
// }
- if (cit->boundaryDualVertex())
- {
- if (cit->featurePointDualVertex())
- {
- boundaryPts[cit->cellIndex()] = featurePoint;
- }
- else if (cit->featureEdgeDualVertex())
- {
- boundaryPts[cit->cellIndex()] = featureEdge;
- }
- else
- {
- 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;
- }
+ boundaryPts[cit->cellIndex()] = classifyBoundaryPoint(cit);
}
else
{
@@ -2040,7 +1499,8 @@ void Foam::conformalVoronoiMesh::indexDualVertices
void Foam::conformalVoronoiMesh::reindexDualVertices
(
- const Map<label>& dualPtIndexMap
+ const Map<label>& dualPtIndexMap,
+ labelList& boundaryPts
)
{
for
@@ -2053,6 +1513,12 @@ void Foam::conformalVoronoiMesh::reindexDualVertices
if (dualPtIndexMap.found(cit->cellIndex()))
{
cit->cellIndex() = dualPtIndexMap[cit->cellIndex()];
+ boundaryPts[cit->cellIndex()] =
+ max
+ (
+ boundaryPts[cit->cellIndex()],
+ boundaryPts[dualPtIndexMap[cit->cellIndex()]]
+ );
}
}
}
@@ -2267,11 +1733,7 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
bool includeEmptyPatches
) const
{
- const label defaultPatchIndex = createPatchInfo
- (
- patchNames,
- patchDicts
- );
+ const label defaultPatchIndex = createPatchInfo(patchNames, patchDicts);
const label nPatches = patchNames.size();
@@ -2296,6 +1758,7 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
List<DynamicList<bool> > indirectPatchFace(nPatches, DynamicList<bool>(0));
+
faces.setSize(number_of_finite_edges());
owner.setSize(number_of_finite_edges());
neighbour.setSize(number_of_finite_edges());
@@ -2765,7 +2228,12 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
// If the two vertices are a pair, then the patch face is
// a desired one.
- if (!isPointPair(vA, vB))
+ if
+ (
+ vA->boundaryPoint() && vB->boundaryPoint()
+ && !ptPairs_.isPointPair(vA, vB)
+ && !ftPtConformer_.featurePointPairs().isPointPair(vA, vB)
+ )
{
indirectPatchFace[patchIndex].append(true);
}
@@ -2797,6 +2265,18 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
if (patchIndex != -1)
{
+// if
+// (
+// vA->boundaryPoint() && vB->boundaryPoint()
+// && !ptPairs_.isPointPair(vA, vB)
+// )
+// {
+// indirectPatchFace[patchIndex].append(true);
+// }
+// else
+// {
+// indirectPatchFace[patchIndex].append(false);
+// }
// patchFaces[patchIndex].append(newDualFace);
// patchOwners[patchIndex].append(own);
// indirectPatchFace[patchIndex].append(false);
@@ -2813,8 +2293,6 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
// {
// patchPPSlaves[patchIndex].append(vA->index());
// }
-
-// baffleFaces[dualFaceI] = patchIndex;
}
// else
{
@@ -2989,7 +2467,6 @@ void Foam::conformalVoronoiMesh::sortProcPatches
faceList& faces = patchFaces[patchI];
labelList& owner = patchOwners[patchI];
DynamicList<label>& slaves = patchPointPairSlaves[patchI];
-
DynamicList<Pair<labelPair> >& sortingIndices
= patchSortingIndices[patchI];
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
index f5a023592c19f373cabf58770d687aa7764bddde..e2b6c0d2e8ebd612a39dcea4ee44f1ca4c16c47f 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C
@@ -62,11 +62,13 @@ void Foam::conformalVoronoiMesh::conformToSurface()
if (Pstream::parRun())
{
- sync(decomposition_().procBounds());
+ sync(decomposition().procBounds());
}
}
else
{
+ ptPairs_.clear();
+
// Rebuild, insert and store new surface conformation
buildSurfaceConformation();
@@ -74,7 +76,7 @@ void Foam::conformalVoronoiMesh::conformToSurface()
{
if (Pstream::parRun())
{
- sync(decomposition_().procBounds());
+ sync(decomposition().procBounds());
}
}
@@ -358,10 +360,13 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
synchroniseSurfaceTrees(surfaceToTreeShape, surfaceHits);
}
+ DynamicList<Vb> pts(2*surfaceHits.size() + 3*featureEdgeHits.size());
+
insertSurfacePointPairs
(
surfaceHits,
- "surfaceConformationLocations_initial.obj"
+ "surfaceConformationLocations_initial.obj",
+ pts
);
// In parallel, synchronise the edge trees
@@ -373,9 +378,21 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
insertEdgePointGroups
(
featureEdgeHits,
- "edgeConformationLocations_initial.obj"
+ "edgeConformationLocations_initial.obj",
+ pts
);
+ pts.shrink();
+
+ Map<label> oldToNewIndices = insertPointPairs(pts, true, true);
+
+ // Re-index the point pairs
+ ptPairs_.reIndex(oldToNewIndices);
+
+ //writePointPairs("pointPairs_initial.obj");
+
+ // Remove location from surface/edge tree
+
timeCheck("After initial conformation");
initialTotalHits = nSurfHits + nFeatEdHits;
@@ -550,215 +567,6 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
}
}
-// for
-// (
-// Delaunay::Finite_edges_iterator eit = finite_edges_begin();
-// eit != finite_edges_end();
-// ++eit
-// )
-// {
-// Cell_handle c = eit->first;
-// Vertex_handle vA = c->vertex(eit->second);
-// Vertex_handle vB = c->vertex(eit->third);
-//
-// if
-// (
-// vA->referred()
-// || vB->referred()
-// )
-// {
-// continue;
-// }
-//
-// if
-// (
-// (vA->internalPoint() && vA->referred())
-// || (vB->internalPoint() && vB->referred())
-// )
-// {
-// continue;
-// }
-//
-// if
-// (
-// (vA->internalPoint() && vB->externalBoundaryPoint())
-// || (vB->internalPoint() && vA->externalBoundaryPoint())
-// || (vA->internalBoundaryPoint() && vB->internalBoundaryPoint())
-// )
-// {
-// pointIndexHitAndFeatureDynList surfaceIntersections(0.5*AtoV);
-// pointIndexHit surfHit;
-// label hitSurface;
-//
-// geometryToConformTo_.findSurfaceNearest
-// (
-// (
-// vA->internalPoint()
-// ? topoint(vA->point())
-// : topoint(vB->point())
-// ),
-// magSqr(topoint(vA->point()) - topoint(vB->point())),
-// surfHit,
-// hitSurface
-// );
-//
-// if (surfHit.hit())
-// {
-// surfaceIntersections.append
-// (
-// pointIndexHitAndFeature(surfHit, hitSurface)
-// );
-//
-// addSurfaceAndEdgeHits
-// (
-// (
-// vA->internalPoint()
-// ? topoint(vA->point())
-// : topoint(vB->point())
-// ),
-// surfaceIntersections,
-// surfacePtReplaceDistCoeffSqr,
-// edgeSearchDistCoeffSqr,
-// surfaceHits,
-// featureEdgeHits,
-// surfaceToTreeShape,
-// edgeToTreeShape,
-// surfacePtToEdgePtDist,
-// false
-// );
-// }
-// }
-// }
-
- for
- (
- Delaunay::Finite_cells_iterator cit = finite_cells_begin();
- cit != finite_cells_end();
- ++cit
- )
- {
- label nInternal = 0;
- for (label vI = 0; vI < 4; vI++)
- {
- if (cit->vertex(vI)->internalPoint())
- {
- nInternal++;
- }
- }
-
- if (nInternal == 1 && cit->hasBoundaryPoint())
- //if (cit->boundaryDualVertex() && !cit->hasReferredPoint())
- {
- const Foam::point& pt = cit->dual();
-
- const Foam::point cellCentre =
- topoint
- (
- CGAL::centroid
- (
- CGAL::Tetrahedron_3<baseK>
- (
- cit->vertex(0)->point(),
- cit->vertex(1)->point(),
- cit->vertex(2)->point(),
- cit->vertex(3)->point()
- )
- )
- );
-
- pointIndexHitAndFeatureDynList surfaceIntersections(0.5*AtoV);
- pointIndexHit surfHit;
- label hitSurface;
-
- geometryToConformTo_.findSurfaceNearestIntersection
- (
- cellCentre,
- pt,
- surfHit,
- hitSurface
- );
-
- if (surfHit.hit())
- {
- surfaceIntersections.append
- (
- pointIndexHitAndFeature(surfHit, hitSurface)
- );
-
- addSurfaceAndEdgeHits
- (
- pt,
- surfaceIntersections,
- surfacePtReplaceDistCoeffSqr,
- edgeSearchDistCoeffSqr,
- surfaceHits,
- featureEdgeHits,
- surfaceToTreeShape,
- edgeToTreeShape,
- surfacePtToEdgePtDist,
- false
- );
- }
- }
- }
-
-// for
-// (
-// Delaunay::Finite_cells_iterator cit = finite_cells_begin();
-// cit != finite_cells_end();
-// ++cit
-// )
-// {
-// if (cit->boundaryDualVertex() && !cit->hasReferredPoint())
-// {
-// const Foam::point& pt = cit->dual();
-//
-// pointIndexHitAndFeatureDynList surfaceIntersections(0.5*AtoV);
-// pointIndexHit surfHit;
-// label hitSurface;
-//
-// geometryToConformTo_.findSurfaceNearest
-// (
-// pt,
-// sqr(0.1*targetCellSize(pt)),
-// surfHit,
-// hitSurface
-// );
-//
-// if (!surfHit.hit())
-// {
-// geometryToConformTo_.findSurfaceNearest
-// (
-// pt,
-// GREAT,
-// surfHit,
-// hitSurface
-// );
-//
-// if (surfHit.hit())
-// {
-// surfaceIntersections.append
-// (
-// pointIndexHitAndFeature(surfHit, hitSurface)
-// );
-//
-// addSurfaceAndEdgeHits
-// (
-// pt,
-// surfaceIntersections,
-// surfacePtReplaceDistCoeffSqr,
-// edgeSearchDistCoeffSqr,
-// surfaceHits,
-// featureEdgeHits,
-// surfaceToTreeShape,
-// edgeToTreeShape,
-// false
-// );
-// }
-// }
-// }
-// }
-
label nVerts = number_of_vertices();
label nSurfHits = surfaceHits.size();
label nFeatEdHits = featureEdgeHits.size();
@@ -789,10 +597,16 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
synchroniseSurfaceTrees(surfaceToTreeShape, surfaceHits);
}
+ DynamicList<Vb> pts
+ (
+ 2*surfaceHits.size() + 3*featureEdgeHits.size()
+ );
+
insertSurfacePointPairs
(
surfaceHits,
- "surfaceConformationLocations_" + name(iterationNo) + ".obj"
+ "surfaceConformationLocations_" + name(iterationNo) + ".obj",
+ pts
);
// In parallel, synchronise the edge trees
@@ -805,9 +619,19 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
insertEdgePointGroups
(
featureEdgeHits,
- "edgeConformationLocations_" + name(iterationNo) + ".obj"
+ "edgeConformationLocations_" + name(iterationNo) + ".obj",
+ pts
);
+ pts.shrink();
+
+ Map<label> oldToNewIndices = insertPointPairs(pts, true, true);
+
+ // Reindex the point pairs
+ ptPairs_.reIndex(oldToNewIndices);
+
+ //writePointPairs("pointPairs_" + name(iterationNo) + ".obj");
+
if (Pstream::parRun())
{
sync
@@ -1155,8 +979,6 @@ bool Foam::conformalVoronoiMesh::dualCellSurfaceAllIntersections
|| is_infinite(fit->first->neighbor(fit->second))
|| !fit->first->hasInternalPoint()
|| !fit->first->neighbor(fit->second)->hasInternalPoint()
-// || fit->first->hasFarPoint()
-// || fit->first->neighbor(fit->second)->hasFarPoint()
)
{
continue;
@@ -1209,14 +1031,13 @@ bool Foam::conformalVoronoiMesh::dualCellSurfaceAllIntersections
const scalar d = p.normalIntersect(r);
- const Foam::point newPoint = vertex + d*n;
+ Foam::point newPoint = vertex + d*n;
pointIndexHitAndFeature info;
-
geometryToConformTo_.findSurfaceNearest
(
newPoint,
- surfaceSearchDistanceSqr(newPoint),
+ 4.0*magSqr(newPoint - vertex),
info.first(),
info.second()
);
@@ -1359,7 +1180,10 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceProtrusion
if
(
is_infinite(c1) || is_infinite(c2)
- || (!c1->hasInternalPoint() && !c2->hasInternalPoint())
+ || (
+ !c1->internalOrBoundaryDualVertex()
+ || !c2->internalOrBoundaryDualVertex()
+ )
|| !c1->real() || !c2->real()
)
{
@@ -1369,19 +1193,24 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceProtrusion
// Foam::point endPt = 0.5*(c1->dual() + c2->dual());
Foam::point endPt = c1->dual();
+ if (magSqr(vert - c1->dual()) < magSqr(vert - c2->dual()))
+ {
+ endPt = c2->dual();
+ }
+
if
(
- magSqr(vert - c1->dual())
- < magSqr(vert - c2->dual())
+ magSqr(vert - endPt)
+ > magSqr(geometryToConformTo().globalBounds().mag())
)
{
- endPt = c2->dual();
+ continue;
}
pointIndexHit surfHit;
label hitSurface;
- geometryToConformTo_.findSurfaceAnyIntersection
+ geometryToConformTo_.findSurfaceNearestIntersection
(
vert,
endPt,
@@ -1401,15 +1230,46 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceProtrusion
const vector& n = norm[0];
- const scalar normalProtrusionDistance =
- (endPt - surfHit.hitPoint()) & n;
+ const scalar normalProtrusionDistance
+ (
+ (endPt - surfHit.hitPoint()) & n
+ );
if (normalProtrusionDistance > maxProtrusionDistance)
{
- surfHitLargest = surfHit;
- hitSurfaceLargest = hitSurface;
+ const plane p(surfHit.hitPoint(), n);
+
+ const plane::ray r(endPt, -n);
+
+ const scalar d = p.normalIntersect(r);
+
+ Foam::point newPoint = endPt - d*n;
+
+ pointIndexHitAndFeature info;
+ geometryToConformTo_.findSurfaceNearest
+ (
+ newPoint,
+ 4.0*magSqr(newPoint - endPt),
+ info.first(),
+ info.second()
+ );
+
+ if (info.first().hit())
+ {
+ if
+ (
+ surfaceLocationConformsToInside
+ (
+ pointIndexHitAndFeature(info.first(), info.second())
+ )
+ )
+ {
+ surfHitLargest = info.first();
+ hitSurfaceLargest = info.second();
- maxProtrusionDistance = normalProtrusionDistance;
+ maxProtrusionDistance = normalProtrusionDistance;
+ }
+ }
}
}
}
@@ -1465,7 +1325,10 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceIncursion
if
(
is_infinite(c1) || is_infinite(c2)
- || (!c1->hasInternalPoint() && !c2->hasInternalPoint())
+ || (
+ !c1->internalOrBoundaryDualVertex()
+ || !c2->internalOrBoundaryDualVertex()
+ )
|| !c1->real() || !c2->real()
)
{
@@ -1475,19 +1338,24 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceIncursion
// Foam::point endPt = 0.5*(c1->dual() + c2->dual());
Foam::point endPt = c1->dual();
+ if (magSqr(vert - c1->dual()) < magSqr(vert - c2->dual()))
+ {
+ endPt = c2->dual();
+ }
+
if
(
- magSqr(vert - c1->dual())
- < magSqr(vert - c2->dual())
+ magSqr(vert - endPt)
+ > magSqr(geometryToConformTo().globalBounds().mag())
)
{
- endPt = c2->dual();
+ continue;
}
pointIndexHit surfHit;
label hitSurface;
- geometryToConformTo_.findSurfaceAnyIntersection
+ geometryToConformTo_.findSurfaceNearestIntersection
(
vert,
endPt,
@@ -1507,19 +1375,46 @@ void Foam::conformalVoronoiMesh::dualCellLargestSurfaceIncursion
const vector& n = norm[0];
- scalar normalIncursionDistance = (endPt - surfHit.hitPoint()) & n;
+ scalar normalIncursionDistance
+ (
+ (endPt - surfHit.hitPoint()) & n
+ );
if (normalIncursionDistance < minIncursionDistance)
{
- surfHitLargest = surfHit;
- hitSurfaceLargest = hitSurface;
+ const plane p(surfHit.hitPoint(), n);
+
+ const plane::ray r(endPt, n);
- minIncursionDistance = normalIncursionDistance;
+ const scalar d = p.normalIntersect(r);
+
+ Foam::point newPoint = endPt + d*n;
+
+ pointIndexHitAndFeature info;
+ geometryToConformTo_.findSurfaceNearest
+ (
+ newPoint,
+ 4.0*magSqr(newPoint - endPt),
+ info.first(),
+ info.second()
+ );
+
+ if (info.first().hit())
+ {
+ if
+ (
+ surfaceLocationConformsToInside
+ (
+ pointIndexHitAndFeature(info.first(), info.second())
+ )
+ )
+ {
+ surfHitLargest = info.first();
+ hitSurfaceLargest = info.second();
- // Info<< nl << "# Incursion: " << endl;
- // meshTools::writeOBJ(Info, vert);
- // meshTools::writeOBJ(Info, edgeMid);
- // Info<< "l Na Nb" << endl;
+ minIncursionDistance = normalIncursionDistance;
+ }
+ }
}
}
}
@@ -2162,9 +2057,11 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
bool isNearFeaturePt = nearFeaturePt(surfPt);
+ bool isNearFeatureEdge = surfacePtNearFeatureEdge(surfPt);
+
bool isNearSurfacePt = nearSurfacePoint(surfHitI);
- if (isNearFeaturePt || isNearSurfacePt)
+ if (isNearFeaturePt || isNearSurfacePt || isNearFeatureEdge)
{
keepSurfacePoint = false;
}
@@ -2219,6 +2116,9 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
// feature edge, give control to edge control points
// instead, this will prevent "pits" forming.
+ // Allow if different surfaces
+
+
keepSurfacePoint = false;
}
@@ -2298,11 +2198,13 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
surfaceHits.append(surfHitI);
appendToSurfacePtTree(surfPt);
surfaceToTreeShape.append(existingSurfacePtLocations_.size() - 1);
+
+// addedPoints.write(surfPt);
+ }
+ else
+ {
+// removedPoints.write(surfPt);
}
-// else
-// {
-// surfaceToTreeShape.remove();
-// }
}
}
@@ -2338,7 +2240,9 @@ void Foam::conformalVoronoiMesh::storeSurfaceConformation()
Vb
(
vit->point(),
- vit->type()
+ vit->index(),
+ vit->type(),
+ Pstream::myProcNo()
)
);
}
@@ -2362,24 +2266,40 @@ void Foam::conformalVoronoiMesh::reinsertSurfaceConformation()
{
Info<< nl << "Reinserting stored surface conformation" << endl;
- const label preReinsertionSize(number_of_vertices());
+ Map<label> oldToNewIndices =
+ insertPointPairs(surfaceConformationVertices_, true, true);
- // It is assumed that the stored surface conformation is on the correct
- // processor and does not need distributed
- rangeInsertWithInfo
- (
- surfaceConformationVertices_.begin(),
- surfaceConformationVertices_.end(),
- true
- );
+ ptPairs_.reIndex(oldToNewIndices);
- const label nInserted = label(number_of_vertices()) - preReinsertionSize;
- const label nFailed = surfaceConformationVertices_.size() - nInserted;
+ PackedBoolList selectedElems(surfaceConformationVertices_.size(), true);
- Info<< " " << returnReduce(nInserted, sumOp<label>())
- << " points reinserted, failed to insert "
- << returnReduce(nFailed, sumOp<label>())
- << endl;
+ forAll(surfaceConformationVertices_, vI)
+ {
+ Vb& v = surfaceConformationVertices_[vI];
+ label& vIndex = v.index();
+
+ Map<label>::const_iterator iter = oldToNewIndices.find(vIndex);
+
+ if (iter != oldToNewIndices.end())
+ {
+ const label newIndex = iter();
+
+ if (newIndex != -1)
+ {
+ vIndex = newIndex;
+ }
+ else
+ {
+ selectedElems[vI] = false;
+ }
+ }
+ }
+
+ inplaceSubset<PackedBoolList, List<Vb> >
+ (
+ selectedElems,
+ surfaceConformationVertices_
+ );
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
index 980064fb48c6e01047a08d9e38bafecd934d3d67..f142470106b2df432714892cc6985ea36afd9bcd 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C
@@ -29,6 +29,7 @@ License
#include "tetrahedron.H"
#include "const_circulator.H"
#include "DelaunayMeshTools.H"
+#include "OBJstream.H"
using namespace Foam::vectorTools;
@@ -225,21 +226,48 @@ void Foam::conformalVoronoiMesh::createEdgePointGroupByCirculating
vector masterPtVec(normalDir + nextNormalDir);
masterPtVec /= mag(masterPtVec) + SMALL;
- if (((normalDir ^ nextNormalDir) & edDir) < SMALL)
+ if
+ (
+ ((normalDir ^ nextNormalDir) & edDir) < SMALL
+ || mag(masterPtVec) < SMALL
+ )
{
// Info<< " IGNORE REGION" << endl;
addedMasterPreviously = false;
- if (circ.size() == 2 && mag((normal & nextNormal) - 1) < SMALL)
+ if
+ (
+ circ.size() == 2
+ && mag((normal & nextNormal) - 1) < SMALL
+ )
{
- // Add an extra point
-
const vector n = 0.5*(normal + nextNormal);
- vector s = ppDist*(edDir ^ n);
-
- createPointPair(ppDist, edgePt + s, n, pts);
- createPointPair(ppDist, edgePt - s, n, pts);
+ const vector s = ppDist*(edDir ^ n);
+
+ if (volType == extendedFeatureEdgeMesh::BOTH)
+ {
+ createBafflePointPair(ppDist, edgePt + s, n, true, pts);
+ createBafflePointPair(ppDist, edgePt - s, n, true, pts);
+ }
+ else
+ {
+ WarningIn
+ (
+ "Foam::conformalVoronoiMesh::"
+ "createEdgePointGroupByCirculating"
+ "("
+ " const extendedFeatureEdgeMesh&,"
+ " const pointIndexHit&,"
+ " DynamicList<Vb>&"
+ ")"
+ ) << "Failed to insert flat/open edge as volType is "
+ << extendedFeatureEdgeMesh::sideVolumeTypeNames_
+ [
+ volType
+ ]
+ << endl;
+ }
break;
}
@@ -247,34 +275,6 @@ void Foam::conformalVoronoiMesh::createEdgePointGroupByCirculating
continue;
}
- if (mag(masterPtVec) < SMALL)
- {
- if (circ.size() == 2)
- {
- // Add an extra point
-
- // Average normal to remove any bias to one side, although as
- // it is a flat edge, the normals should be essentially the same
- const vector n = 0.5*(normal + nextNormal);
-
- // Direction along the surface to the control point, sense of
- // direction not important, as +s and -s can be used because it
- // is a flat edge
- vector s = ppDist*(edDir ^ n);
-
- createPointPair(ppDist, edgePt + s, n, pts);
- createPointPair(ppDist, edgePt - s, n, pts);
-
- break;
- }
- else
- {
-// Info<< " IGNORE REGION" << endl;
- addedMasterPreviously = false;
- continue;
- }
- }
-
const Foam::point masterPt = edgePt + ppDist*masterPtVec;
// Check that region is inside or outside
@@ -490,11 +490,19 @@ void Foam::conformalVoronoiMesh::createExternalEdgePointGroup
const vectorField& feNormals = feMesh.normals();
const labelList& edNormalIs = feMesh.edgeNormals()[edHit.index()];
+ const List<extendedFeatureEdgeMesh::sideVolumeType>& normalVolumeTypes =
+ feMesh.normalVolumeTypes();
// As this is an external edge, there are two normals by definition
const vector& nA = feNormals[edNormalIs[0]];
const vector& nB = feNormals[edNormalIs[1]];
+ const extendedFeatureEdgeMesh::sideVolumeType& volTypeA =
+ normalVolumeTypes[edNormalIs[0]];
+
+ const extendedFeatureEdgeMesh::sideVolumeType& volTypeB =
+ normalVolumeTypes[edNormalIs[1]];
+
if (areParallel(nA, nB))
{
// The normals are nearly parallel, so this is too sharp a feature to
@@ -524,13 +532,6 @@ void Foam::conformalVoronoiMesh::createExternalEdgePointGroup
// Convex. So refPt will be inside domain and hence a master point
Foam::point refPt = edgePt - ppDist*refVec;
- // Result when the points are eventually inserted.
- // Add number_of_vertices() at insertion of first vertex to all numbers:
- // pt index type
- // refPt 0 1
- // reflectedA 1 0
- // reflectedB 2 0
-
// Insert the master point pairing the the first slave
if (!geometryToConformTo_.inside(refPt))
@@ -559,7 +560,11 @@ void Foam::conformalVoronoiMesh::createExternalEdgePointGroup
(
reflectedA,
vertexCount() + pts.size(),
- Vb::vtExternalFeatureEdge,
+ (
+ volTypeA == extendedFeatureEdgeMesh::BOTH
+ ? Vb::vtInternalFeatureEdge
+ : Vb::vtExternalFeatureEdge
+ ),
Pstream::myProcNo()
)
);
@@ -571,10 +576,26 @@ void Foam::conformalVoronoiMesh::createExternalEdgePointGroup
(
reflectedB,
vertexCount() + pts.size(),
- Vb::vtExternalFeatureEdge,
+ (
+ volTypeB == extendedFeatureEdgeMesh::BOTH
+ ? Vb::vtInternalFeatureEdge
+ : Vb::vtExternalFeatureEdge
+ ),
Pstream::myProcNo()
)
);
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(),
+ pts[pts.size() - 1].index()
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(),
+ pts[pts.size() - 2].index()
+ );
}
@@ -591,11 +612,19 @@ void Foam::conformalVoronoiMesh::createInternalEdgePointGroup
const vectorField& feNormals = feMesh.normals();
const labelList& edNormalIs = feMesh.edgeNormals()[edHit.index()];
+ const List<extendedFeatureEdgeMesh::sideVolumeType>& normalVolumeTypes =
+ feMesh.normalVolumeTypes();
// As this is an external edge, there are two normals by definition
const vector& nA = feNormals[edNormalIs[0]];
const vector& nB = feNormals[edNormalIs[1]];
+ const extendedFeatureEdgeMesh::sideVolumeType& volTypeA =
+ normalVolumeTypes[edNormalIs[0]];
+
+// const extendedFeatureEdgeMesh::sideVolumeType& volTypeB =
+// normalVolumeTypes[edNormalIs[1]];
+
if (areParallel(nA, nB))
{
// The normals are nearly parallel, so this is too sharp a feature to
@@ -705,14 +734,30 @@ void Foam::conformalVoronoiMesh::createInternalEdgePointGroup
(
reflMasterPt,
vertexCount() + pts.size(),
- Vb::vtExternalFeatureEdge,
+ (
+ volTypeA == extendedFeatureEdgeMesh::BOTH
+ ? Vb::vtInternalFeatureEdge
+ : Vb::vtExternalFeatureEdge
+ ),
Pstream::myProcNo()
)
);
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(),
+ pts[pts.size() - 1].index()
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(),
+ pts[pts.size() - 1].index()
+ );
+
if (nAddPoints == 1)
{
- // One additinal point is the reflection of the slave point,
+ // One additional point is the reflection of the slave point,
// i.e. the original reference point
pts.append
(
@@ -767,6 +812,8 @@ void Foam::conformalVoronoiMesh::createFlatEdgePointGroup
const vectorField& feNormals = feMesh.normals();
const labelList& edNormalIs = feMesh.edgeNormals()[edHit.index()];
+ const List<extendedFeatureEdgeMesh::sideVolumeType>& normalVolumeTypes =
+ feMesh.normalVolumeTypes();
// As this is a flat edge, there are two normals by definition
const vector& nA = feNormals[edNormalIs[0]];
@@ -781,8 +828,21 @@ void Foam::conformalVoronoiMesh::createFlatEdgePointGroup
// is a flat edge
vector s = ppDist*(feMesh.edgeDirections()[edHit.index()] ^ n);
- createPointPair(ppDist, edgePt + s, n, pts);
- createPointPair(ppDist, edgePt - s, n, pts);
+ if (normalVolumeTypes[edNormalIs[0]] == extendedFeatureEdgeMesh::OUTSIDE)
+ {
+ createPointPair(ppDist, edgePt + s, -n, true, pts);
+ createPointPair(ppDist, edgePt - s, -n, true, pts);
+ }
+ else if (normalVolumeTypes[edNormalIs[0]] == extendedFeatureEdgeMesh::BOTH)
+ {
+ createBafflePointPair(ppDist, edgePt + s, n, true, pts);
+ createBafflePointPair(ppDist, edgePt - s, n, true, pts);
+ }
+ else
+ {
+ createPointPair(ppDist, edgePt + s, n, true, pts);
+ createPointPair(ppDist, edgePt - s, n, true, pts);
+ }
}
@@ -814,16 +874,23 @@ void Foam::conformalVoronoiMesh::createOpenEdgePointGroup
plane facePlane(edgePt, n);
- Foam::point pt1 = edgePt + s + ppDist*n;
- Foam::point pt2 = edgePt - s + ppDist*n;
+ const label initialPtsSize = pts.size();
- Foam::point pt3 = facePlane.mirror(pt1);
- Foam::point pt4 = facePlane.mirror(pt2);
+ if
+ (
+ !geometryToConformTo_.inside(edgePt)
+ )
+ {
+ return;
+ }
+
+ createBafflePointPair(ppDist, edgePt - s, n, true, pts);
+ createBafflePointPair(ppDist, edgePt + s, n, false, pts);
- pts.append(Vb(pt1, Vb::vtInternalSurface));
- pts.append(Vb(pt2, Vb::vtInternalSurface));
- pts.append(Vb(pt3, Vb::vtInternalSurface));
- pts.append(Vb(pt4, Vb::vtInternalSurface));
+ for (label ptI = initialPtsSize; ptI < pts.size(); ++ptI)
+ {
+ pts[ptI].type() = Vb::vtInternalFeatureEdge;
+ }
}
else
{
@@ -846,26 +913,245 @@ void Foam::conformalVoronoiMesh::createMultipleEdgePointGroup
const scalar ppDist = pointPairDistance(edgePt);
+ const vector edDir = feMesh.edgeDirections()[edHit.index()];
+
const vectorField& feNormals = feMesh.normals();
const labelList& edNormalIs = feMesh.edgeNormals()[edHit.index()];
+ const labelList& normalDirs = feMesh.normalDirections()[edHit.index()];
- Info<< edNormalIs.size() << endl;
+ const List<extendedFeatureEdgeMesh::sideVolumeType>& normalVolumeTypes =
+ feMesh.normalVolumeTypes();
- // As this is a flat edge, there are two normals by definition
- const vector& nA = feNormals[edNormalIs[0]];
- const vector& nB = feNormals[edNormalIs[1]];
+ labelList nNormalTypes(4, 0);
- // Average normal to remove any bias to one side, although as this
- // is a flat edge, the normals should be essentially the same
- const vector n = 0.5*(nA + nB);
+ forAll(edNormalIs, edgeNormalI)
+ {
+ const extendedFeatureEdgeMesh::sideVolumeType sType =
+ normalVolumeTypes[edNormalIs[edgeNormalI]];
- // Direction along the surface to the control point, sense of edge
- // direction not important, as +s and -s can be used because this
- // is a flat edge
- vector s = ppDist*(feMesh.edgeDirections()[edHit.index()] ^ n);
+ nNormalTypes[sType]++;
+ }
+
+ if (nNormalTypes[extendedFeatureEdgeMesh::BOTH] == 4)
+ {
+ label masterEdgeNormalIndex = -1;
- createPointPair(ppDist, edgePt + s, n, pts);
- createPointPair(ppDist, edgePt - s, n, pts);
+ forAll(edNormalIs, edgeNormalI)
+ {
+ const extendedFeatureEdgeMesh::sideVolumeType sType =
+ normalVolumeTypes[edNormalIs[edgeNormalI]];
+
+ if (sType == extendedFeatureEdgeMesh::BOTH)
+ {
+ masterEdgeNormalIndex = edgeNormalI;
+ break;
+ }
+ }
+
+ const vector& n = feNormals[edNormalIs[masterEdgeNormalIndex]];
+
+ label nDir = normalDirs[masterEdgeNormalIndex];
+
+ vector normalDir =
+ (feNormals[edNormalIs[masterEdgeNormalIndex]] ^ edDir);
+ normalDir *= nDir/mag(normalDir);
+
+ Foam::point pt1 = edgePt + ppDist*normalDir + ppDist*n;
+ Foam::point pt2 = edgePt + ppDist*normalDir - ppDist*n;
+
+ plane plane3(edgePt, normalDir);
+
+ Foam::point pt3 = plane3.mirror(pt1);
+ Foam::point pt4 = plane3.mirror(pt2);
+
+ pts.append
+ (
+ Vb
+ (
+ pt1,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+ pts.append
+ (
+ Vb
+ (
+ pt2,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+
+ pts.append
+ (
+ Vb
+ (
+ pt3,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+
+ pts.append
+ (
+ Vb
+ (
+ pt4,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+ }
+ else if
+ (
+ nNormalTypes[extendedFeatureEdgeMesh::BOTH] == 1
+ && nNormalTypes[extendedFeatureEdgeMesh::INSIDE] == 2
+ )
+ {
+ label masterEdgeNormalIndex = -1;
+
+ forAll(edNormalIs, edgeNormalI)
+ {
+ const extendedFeatureEdgeMesh::sideVolumeType sType =
+ normalVolumeTypes[edNormalIs[edgeNormalI]];
+
+ if (sType == extendedFeatureEdgeMesh::BOTH)
+ {
+ masterEdgeNormalIndex = edgeNormalI;
+ break;
+ }
+ }
+
+ const vector& n = feNormals[edNormalIs[masterEdgeNormalIndex]];
+
+ label nDir = normalDirs[masterEdgeNormalIndex];
+
+ vector normalDir =
+ (feNormals[edNormalIs[masterEdgeNormalIndex]] ^ edDir);
+ normalDir *= nDir/mag(normalDir);
+
+ const label nextNormalI =
+ (masterEdgeNormalIndex + 1) % edNormalIs.size();
+ if ((normalDir & feNormals[edNormalIs[nextNormalI]]) > 0)
+ {
+ normalDir *= -1;
+ }
+
+ Foam::point pt1 = edgePt + ppDist*normalDir + ppDist*n;
+ Foam::point pt2 = edgePt + ppDist*normalDir - ppDist*n;
+
+ plane plane3(edgePt, normalDir);
+
+ Foam::point pt3 = plane3.mirror(pt1);
+ Foam::point pt4 = plane3.mirror(pt2);
+
+ pts.append
+ (
+ Vb
+ (
+ pt1,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+ pts.append
+ (
+ Vb
+ (
+ pt2,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+
+ pts.append
+ (
+ Vb
+ (
+ pt3,
+ vertexCount() + pts.size(),
+ Vb::vtExternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+
+ pts.append
+ (
+ Vb
+ (
+ pt4,
+ vertexCount() + pts.size(),
+ Vb::vtExternalSurface,
+ Pstream::myProcNo()
+ )
+ );
+
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 3].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
+ }
+
+
+// // As this is a flat edge, there are two normals by definition
+// const vector& nA = feNormals[edNormalIs[0]];
+// const vector& nB = feNormals[edNormalIs[1]];
+//
+// // Average normal to remove any bias to one side, although as this
+// // is a flat edge, the normals should be essentially the same
+// const vector n = 0.5*(nA + nB);
+//
+// // Direction along the surface to the control point, sense of edge
+// // direction not important, as +s and -s can be used because this
+// // is a flat edge
+// vector s = ppDist*(feMesh.edgeDirections()[edHit.index()] ^ n);
+//
+// createBafflePointPair(ppDist, edgePt + s, n, true, pts);
+// createBafflePointPair(ppDist, edgePt - s, n, true, pts);
}
@@ -884,7 +1170,10 @@ void Foam::conformalVoronoiMesh::insertFeaturePoints(bool distribute)
const List<Vb>& ftPtVertices = ftPtConformer_.featurePointVertices();
// Insert the created points directly as already distributed.
- this->DelaunayMesh<Delaunay>::insertPoints(ftPtVertices);
+ Map<label> oldToNewIndices =
+ this->DelaunayMesh<Delaunay>::insertPoints(ftPtVertices, true);
+
+ ftPtConformer_.reIndexPointPairs(oldToNewIndices);
label nFeatureVertices = number_of_vertices() - preFeaturePointSize;
reduce(nFeatureVertices, sumOp<label>());
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
index e9e03a63304153cf8478a8c2f08a9a45bd79542a..7d35594258dead100e872f25c60909fec799b2f4 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshI.H
@@ -225,6 +225,7 @@ inline void Foam::conformalVoronoiMesh::createPointPair
const scalar ppDist,
const Foam::point& surfPt,
const vector& n,
+ const bool ptPair,
DynamicList<Vb>& pts
) const
{
@@ -260,14 +261,14 @@ inline void Foam::conformalVoronoiMesh::createPointPair
)
);
-// if (ptPair)
-// {
-// ptPairs_.addPointPair
-// (
-// pts[pts.size() - 2].index(),
-// pts[pts.size() - 1].index() // external 0 -> slave
-// );
-// }
+ if (ptPair)
+ {
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(),
+ pts[pts.size() - 1].index() // external 0 -> slave
+ );
+ }
}
// else
// {
@@ -305,6 +306,7 @@ inline void Foam::conformalVoronoiMesh::createBafflePointPair
const scalar ppDist,
const Foam::point& surfPt,
const vector& n,
+ const bool ptPair,
DynamicList<Vb>& pts
) const
{
@@ -315,8 +317,8 @@ inline void Foam::conformalVoronoiMesh::createBafflePointPair
Vb
(
surfPt - ppDistn,
- vertexCount() + 1 + pts.size(),
- Vb::vtInternalSurface,
+ vertexCount() + pts.size(),
+ Vb::vtInternalSurfaceBaffle,
Pstream::myProcNo()
)
);
@@ -326,96 +328,35 @@ inline void Foam::conformalVoronoiMesh::createBafflePointPair
Vb
(
surfPt + ppDistn,
- vertexCount() + 1 + pts.size(),
- Vb::vtInternalSurface,
+ vertexCount() + pts.size(),
+ Vb::vtExternalSurfaceBaffle,
Pstream::myProcNo()
)
);
-}
-
-inline bool Foam::conformalVoronoiMesh::internalPointIsInside
-(
- const Foam::point& pt
-) const
-{
- if
- (
- !geometryToConformTo_.inside(pt)
- || !geometryToConformTo_.globalBounds().contains(pt)
- )
+ if (ptPair)
{
- return false;
+ ptPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(), // external 0 -> slave
+ pts[pts.size() - 1].index()
+ );
}
-
- return true;
}
-inline bool Foam::conformalVoronoiMesh::isPointPair
+inline bool Foam::conformalVoronoiMesh::internalPointIsInside
(
- const Vertex_handle& vA,
- const Vertex_handle& vB
+ const Foam::point& pt
) const
{
- // Want to do this topologically, but problem if vertices are redistributed
- // so that one of the point pair is one processor and the other is on
- // another.
-
- const Foam::point& ptA = topoint(vA->point());
- const Foam::point& ptB = topoint(vB->point());
-
if
(
- (
- vA->type() == Vb::vtInternalSurface
- && vB->type() == Vb::vtExternalSurface
- )
- ||
- (
- vB->type() == Vb::vtInternalSurface
- && vA->type() == Vb::vtExternalSurface
- )
- ||
- (
- vA->type() == Vb::vtInternalFeatureEdge
- && vB->type() == Vb::vtExternalFeatureEdge
- )
- ||
- (
- vB->type() == Vb::vtInternalFeatureEdge
- && vA->type() == Vb::vtExternalFeatureEdge
- )
- ||
- (
- vA->type() == Vb::vtInternalSurface
- && vB->type() == Vb::vtExternalFeatureEdge
- )
- ||
- (
- vB->type() == Vb::vtInternalSurface
- && vA->type() == Vb::vtExternalFeatureEdge
- )
- ||
- (
- vA->type() == Vb::vtExternalSurface
- && vB->type() == Vb::vtInternalFeatureEdge
- )
- ||
- (
- vB->type() == Vb::vtExternalSurface
- && vA->type() == Vb::vtInternalFeatureEdge
- )
+ !geometryToConformTo_.globalBounds().contains(pt)
+ || !geometryToConformTo_.inside(pt)
)
{
- const scalar distSqr = magSqr(ptA - ptB);
-
- const scalar ppDistSqr = sqr(2*pointPairDistance(0.5*(ptA + ptB)));
-
- if (distSqr > 1.001*ppDistSqr)
- {
- return false;
- }
+ return false;
}
return true;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C
index d733b33e35ae9c623f1abf898a192f707efae3f8..b0da21cbbb12b2cef3e4a8316d72d81129edb107 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C
@@ -36,6 +36,8 @@ License
#include "indexedVertexOps.H"
#include "DelaunayMeshTools.H"
#include "syncTools.H"
+#include "faceSet.H"
+#include "OBJstream.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@@ -114,7 +116,6 @@ void Foam::conformalVoronoiMesh::writeMesh(const fileName& instance)
wordList patchNames;
PtrList<dictionary> patchDicts;
pointField cellCentres;
-
PackedBoolList boundaryFacesToRemove;
calcDualMesh
@@ -789,7 +790,7 @@ void Foam::conformalVoronoiMesh::writeMesh
const wordList& patchNames,
const PtrList<dictionary>& patchDicts,
const pointField& cellCentres,
- const PackedBoolList& boundaryFacesToRemove
+ PackedBoolList& boundaryFacesToRemove
) const
{
if (foamyHexMeshControls().objOutput())
@@ -912,42 +913,84 @@ void Foam::conformalVoronoiMesh::writeMesh
- // Add indirectPatchFaces to a face zone
+ Info<< indent << "Add pointZones" << endl;
{
- labelList addr(boundaryFacesToRemove.count());
- label count = 0;
+ label sz = mesh.pointZones().size();
- forAll(boundaryFacesToRemove, faceI)
+ DynamicList<label> bPts(boundaryPts.size());
+
+ forAll(dualMeshPointTypeNames_, typeI)
{
- if
- (
- boundaryFacesToRemove[faceI]
- && mesh.faceZones().whichZone(faceI) == -1
- )
+ forAll(boundaryPts, ptI)
{
- addr[count++] = faceI;
+ const label& bPtType = boundaryPts[ptI];
+
+ if (bPtType == typeI)
+ {
+ bPts.append(ptI);
+ }
}
- }
- addr.setSize(count);
+// syncTools::syncPointList(mesh, bPts, maxEqOp<label>(), -1);
- label sz = mesh.faceZones().size();
- boolList flip(addr.size(), false);
- mesh.faceZones().setSize(sz + 1);
- mesh.faceZones().set
- (
- sz,
- new faceZone
+ Info<< incrIndent << indent
+ << "Adding " << bPts.size()
+ << " points of type " << dualMeshPointTypeNames_.words()[typeI]
+ << decrIndent << endl;
+
+ mesh.pointZones().append
(
- "indirectPatchFaces",
- addr,
- flip,
- sz,
- mesh.faceZones()
- )
- );
+ new pointZone
+ (
+ dualMeshPointTypeNames_.words()[typeI],
+ bPts,
+ sz + typeI,
+ mesh.pointZones()
+ )
+ );
+
+ bPts.clear();
+ }
+ }
+
+
+
+ // Add indirectPatchFaces to a face zone
+ Info<< indent << "Adding indirect patch faces set" << endl;
+
+ syncTools::syncFaceList
+ (
+ mesh,
+ boundaryFacesToRemove,
+ orOp<bool>()
+ );
+
+ labelList addr(boundaryFacesToRemove.count());
+ label count = 0;
+
+ forAll(boundaryFacesToRemove, faceI)
+ {
+ if (boundaryFacesToRemove[faceI])
+ {
+ addr[count++] = faceI;
+ }
}
+ addr.setSize(count);
+
+ faceSet indirectPatchFaces
+ (
+ mesh,
+ "indirectPatchFaces",
+ addr,
+ IOobject::AUTO_WRITE
+ );
+
+ indirectPatchFaces.sync(mesh);
+
+
+ Info<< decrIndent;
+
timeCheck("Before fvMesh filtering");
autoPtr<polyMeshFilter> meshFilter;
@@ -966,9 +1009,11 @@ void Foam::conformalVoronoiMesh::writeMesh
{
const autoPtr<fvMesh>& newMesh = meshFilter().filteredMesh();
- polyTopoChange meshMod(newMesh);
+ polyTopoChange meshMod(newMesh());
- meshMod.changeMesh(mesh, false);
+ autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
+
+ polyMeshFilter::copySets(newMesh(), mesh);
}
}
@@ -1001,9 +1046,11 @@ void Foam::conformalVoronoiMesh::writeMesh
{
const autoPtr<fvMesh>& newMesh = meshFilter().filteredMesh();
- polyTopoChange meshMod(newMesh);
+ polyTopoChange meshMod(newMesh());
+
+ autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
- meshMod.changeMesh(mesh, false);
+ polyMeshFilter::copySets(newMesh(), mesh);
}
}
@@ -1068,7 +1115,7 @@ void Foam::conformalVoronoiMesh::writeMesh
// writeCellCentres(mesh);
-// findRemainingProtrusionSet(mesh);
+ findRemainingProtrusionSet(mesh);
}
@@ -1375,4 +1422,33 @@ Foam::labelHashSet Foam::conformalVoronoiMesh::findRemainingProtrusionSet
}
+void Foam::conformalVoronoiMesh::writePointPairs
+(
+ const fileName& fName
+) const
+{
+ OBJstream os(fName);
+
+ for
+ (
+ Delaunay::Finite_edges_iterator eit = finite_edges_begin();
+ eit != finite_edges_end();
+ ++eit
+ )
+ {
+ Cell_handle c = eit->first;
+ Vertex_handle vA = c->vertex(eit->second);
+ Vertex_handle vB = c->vertex(eit->third);
+
+ if (ptPairs_.isPointPair(vA, vB))
+ {
+ os.write
+ (
+ linePointRef(topoint(vA->point()), topoint(vB->point()))
+ );
+ }
+ }
+}
+
+
// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.C
index 542f84afe52657e9bfe87c9cc7e4152729ac5be0..fa42d63e1737644cea9661048d4677c6d3b5166c 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.C
@@ -146,7 +146,7 @@ void Foam::featurePointConformer::addMasterAndSlavePoints
)
);
-// const label masterIndex = pts[pts.size() - 1].index();
+ const label masterIndex = pts.last().index();
//meshTools::writeOBJ(strMasters, masterPt);
@@ -180,6 +180,11 @@ void Foam::featurePointConformer::addMasterAndSlavePoints
)
);
+ ftPtPairs_.addPointPair
+ (
+ masterIndex,
+ pts.last().index()
+ );
//meshTools::writeOBJ(strSlaves, slavePt);
}
@@ -525,7 +530,8 @@ Foam::featurePointConformer::featurePointConformer
foamyHexMesh_(foamyHexMesh),
foamyHexMeshControls_(foamyHexMesh.foamyHexMeshControls()),
geometryToConformTo_(foamyHexMesh.geometryToConformTo()),
- featurePointVertices_()
+ featurePointVertices_(),
+ ftPtPairs_(foamyHexMesh)
{
Info<< nl
<< "Conforming to feature points" << endl;
@@ -598,6 +604,30 @@ void Foam::featurePointConformer::distribute
{
featurePointVertices_[vI].procIndex() = Pstream::myProcNo();
}
+
+ // Also update the feature point pairs
+}
+
+
+void Foam::featurePointConformer::reIndexPointPairs
+(
+ const Map<label>& oldToNewIndices
+)
+{
+ forAll(featurePointVertices_, vI)
+ {
+ const label currentIndex = featurePointVertices_[vI].index();
+
+ Map<label>::const_iterator newIndexIter =
+ oldToNewIndices.find(currentIndex);
+
+ if (newIndexIter != oldToNewIndices.end())
+ {
+ featurePointVertices_[vI].index() = newIndexIter();
+ }
+ }
+
+ ftPtPairs_.reIndex(oldToNewIndices);
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.H
index d4d52c1fff88eea251fef5bc0e368c15b4ce8c9b..e3c415b7fabe3997d1b7ecb8d8b7313b2b262916 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformer.H
@@ -44,6 +44,7 @@ SourceFiles
#include "DynamicList.H"
#include "List.H"
#include "extendedFeatureEdgeMesh.H"
+#include "pointPairs.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -83,6 +84,9 @@ class featurePointConformer
// triangulation clear.
List<Vb> featurePointVertices_;
+ //-
+ mutable pointPairs<Delaunay> ftPtPairs_;
+
// Private Member Functions
@@ -163,12 +167,18 @@ public:
// triangulation.
inline const List<Vb>& featurePointVertices() const;
+ //- Return the feature point pair table
+ inline const pointPairs<Delaunay>& featurePointPairs() const;
+
// Edit
//- Distribute the feature point vertices according to the
// supplied background mesh
void distribute(const backgroundMeshDecomposition& decomposition);
+
+ //- Reindex the feature point pairs using the map.
+ void reIndexPointPairs(const Map<label>& oldToNewIndices);
};
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerI.H
index 406ed3f00048ce0890e612e901f82df573ebd377..0666f274ba0b678d31a149eccd8a47b85d84b584 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerI.H
@@ -28,5 +28,11 @@ const Foam::List<Vb>& Foam::featurePointConformer::featurePointVertices() const
return featurePointVertices_;
}
+const Foam::pointPairs<Delaunay>&
+Foam::featurePointConformer::featurePointPairs() const
+{
+ return ftPtPairs_;
+}
+
// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerSpecialisations.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerSpecialisations.C
index 5ddccc0d68a3cbdcead08d9c961f8d9febd0abe2..b12f56b909e49a829bd5a2c2d13b6d999adca936 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerSpecialisations.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/featurePointConformer/featurePointConformerSpecialisations.C
@@ -168,9 +168,17 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtA, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ internalPtA,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
);
+ const label internalPtAIndex(pts.last().index());
+
const Foam::point internalPtB =
concaveEdgeExternalPt
- 2.0*planeB.distance(concaveEdgeExternalPt)
@@ -178,9 +186,17 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtB, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ internalPtB,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
);
+ const label internalPtBIndex(pts.last().index());
+
// Add the external points
Foam::point externalPtD;
@@ -288,7 +304,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtD, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ externalPtD,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtAIndex,
+ pts.last().index()
);
}
}
@@ -319,7 +347,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtE, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ externalPtE,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtBIndex,
+ pts.last().index()
);
}
}
@@ -328,9 +368,29 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(concaveEdgeExternalPt, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ concaveEdgeExternalPt,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
);
+ ftPtPairs_.addPointPair
+ (
+ internalPtBIndex,
+ pts.last().index()
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtAIndex,
+ pts.last().index()
+ );
+
+ const label concaveEdgeExternalPtIndex(pts.last().index());
+
const scalar totalAngle = radToDeg
(
constant::mathematical::pi
@@ -377,7 +437,21 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtF, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ internalPtF,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ const label internalPtFIndex(pts.last().index());
+
+ ftPtPairs_.addPointPair
+ (
+ concaveEdgeExternalPtIndex,
+ pts.last().index()
);
const Foam::point externalPtG =
@@ -386,7 +460,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtG, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ externalPtG,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtFIndex,
+ pts.last().index()
);
}
@@ -520,9 +606,17 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtA, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ internalPtA,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
);
+ const label internalPtAIndex(pts.last().index());
+
const Foam::point internalPtB =
convexEdgeExternalPt
+ 2.0*planeB.distance(convexEdgeExternalPt)
@@ -530,9 +624,17 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtB, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ internalPtB,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
);
+ const label internalPtBIndex(pts.last().index());
+
// Add the internal points
Foam::point externalPtD;
@@ -643,7 +745,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtD, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ externalPtD,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtAIndex,
+ pts.last().index()
);
}
}
@@ -674,7 +788,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtE, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ externalPtE,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtBIndex,
+ pts.last().index()
);
}
}
@@ -683,7 +809,25 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(convexEdgeExternalPt, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ convexEdgeExternalPt,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtBIndex,
+ pts.last().index()
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ internalPtAIndex,
+ pts.last().index()
);
const scalar totalAngle = radToDeg
@@ -732,7 +876,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(internalPtF, Vb::vtExternalFeaturePoint)
+ Vb
+ (
+ internalPtF,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtExternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(),
+ pts.last().index()
);
const Foam::point externalPtG =
@@ -741,7 +897,19 @@ bool Foam::featurePointConformer::createSpecialisedFeaturePoint
pts.append
(
- Vb(externalPtG, Vb::vtInternalFeaturePoint)
+ Vb
+ (
+ externalPtG,
+ foamyHexMesh_.vertexCount() + pts.size(),
+ Vb::vtInternalFeaturePoint,
+ Pstream::myProcNo()
+ )
+ );
+
+ ftPtPairs_.addPointPair
+ (
+ pts[pts.size() - 2].index(),
+ pts.last().index()
);
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
index 15b155d85b1f92aa131bdd8a9a31c060aed28bb1..f4397070c9a333e59fdb717e4673709502e9442f 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCell.H
@@ -193,6 +193,8 @@ public:
//- Does the Dual vertex form part of a processor patch
inline bool parallelDualVertex() const;
+ inline Foam::label vertexLowestProc() const;
+
//- Using the globalIndex object, return a list of four (sorted) global
// Delaunay vertex indices that uniquely identify this tet in parallel
inline Foam::tetCell vertexGlobalIndices
@@ -219,7 +221,9 @@ public:
// least one Delaunay vertex outside and at least one inside
inline bool boundaryDualVertex() const;
- inline bool baffleBoundaryDualVertex() const;
+ inline bool baffleSurfaceDualVertex() const;
+
+ inline bool baffleEdgeDualVertex() const;
//- A dual vertex on a feature edge will result from this Delaunay cell
inline bool featureEdgeDualVertex() const;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
index 6c55a8eb30706bc938a4cdb02ccefb25600fec02..d0ae6137a5d33211d65581ee86a22b0ee31f0722 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedCell/indexedCellI.H
@@ -291,6 +291,23 @@ inline bool CGAL::indexedCell<Gt, Cb>::parallelDualVertex() const
}
+template<class Gt, class Cb>
+inline Foam::label CGAL::indexedCell<Gt, Cb>::vertexLowestProc() const
+{
+ Foam::label lowestProc = -1;
+
+ for (int i = 0; i < 4; ++i)
+ {
+ if (this->vertex(i)->referred())
+ {
+ lowestProc = min(lowestProc, this->vertex(i)->procIndex());
+ }
+ }
+
+ return lowestProc;
+}
+
+
template<class Gt, class Cb>
inline Foam::tetCell CGAL::indexedCell<Gt, Cb>::vertexGlobalIndices
(
@@ -426,21 +443,42 @@ inline bool CGAL::indexedCell<Gt, Cb>::boundaryDualVertex() const
template<class Gt, class Cb>
-inline bool CGAL::indexedCell<Gt, Cb>::baffleBoundaryDualVertex() const
+inline bool CGAL::indexedCell<Gt, Cb>::baffleSurfaceDualVertex() const
+{
+ return
+ (
+ (
+ this->vertex(0)->internalBaffleSurfacePoint()
+ || this->vertex(1)->internalBaffleSurfacePoint()
+ || this->vertex(2)->internalBaffleSurfacePoint()
+ || this->vertex(3)->internalBaffleSurfacePoint()
+ )
+ && (
+ this->vertex(0)->externalBaffleSurfacePoint()
+ || this->vertex(1)->externalBaffleSurfacePoint()
+ || this->vertex(2)->externalBaffleSurfacePoint()
+ || this->vertex(3)->externalBaffleSurfacePoint()
+ )
+ );
+}
+
+
+template<class Gt, class Cb>
+inline bool CGAL::indexedCell<Gt, Cb>::baffleEdgeDualVertex() const
{
return
(
(
- this->vertex(0)->internalBafflePoint()
- || this->vertex(1)->internalBafflePoint()
- || this->vertex(2)->internalBafflePoint()
- || this->vertex(3)->internalBafflePoint()
+ this->vertex(0)->internalBaffleEdgePoint()
+ || this->vertex(1)->internalBaffleEdgePoint()
+ || this->vertex(2)->internalBaffleEdgePoint()
+ || this->vertex(3)->internalBaffleEdgePoint()
)
&& (
- this->vertex(0)->externalBafflePoint()
- || this->vertex(1)->externalBafflePoint()
- || this->vertex(2)->externalBafflePoint()
- || this->vertex(3)->externalBafflePoint()
+ this->vertex(0)->externalBaffleEdgePoint()
+ || this->vertex(1)->externalBaffleEdgePoint()
+ || this->vertex(2)->externalBaffleEdgePoint()
+ || this->vertex(3)->externalBaffleEdgePoint()
)
);
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.H
index 362e457d74895ce0492a33bbd1b76b499ff9ed0e..311329da0a9db34a870aec4072ff8948a3bf857a 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.H
@@ -243,10 +243,12 @@ public:
inline bool surfacePoint() const;
inline bool internalBoundaryPoint() const;
- inline bool internalBafflePoint() const;
+ inline bool internalBaffleSurfacePoint() const;
+ inline bool internalBaffleEdgePoint() const;
inline bool externalBoundaryPoint() const;
- inline bool externalBafflePoint() const;
+ inline bool externalBaffleSurfacePoint() const;
+ inline bool externalBaffleEdgePoint() const;
inline bool constrained() const;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexI.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexI.H
index c5930395c3df929a497d94cd085020aa1f2cb1f7..688e9f28e1e27f3598e345ceb5ebd6483557349b 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexI.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertexI.H
@@ -308,15 +308,16 @@ inline bool CGAL::indexedVertex<Gt, Vb>::internalBoundaryPoint() const
}
template<class Gt, class Vb>
-inline bool CGAL::indexedVertex<Gt, Vb>::internalBafflePoint() const
+inline bool CGAL::indexedVertex<Gt, Vb>::internalBaffleSurfacePoint() const
{
- return
- (
- type_ == vtInternalSurfaceBaffle
- || type_ == vtInternalFeatureEdgeBaffle
- );
+ return (type_ == vtInternalSurfaceBaffle);
}
+template<class Gt, class Vb>
+inline bool CGAL::indexedVertex<Gt, Vb>::internalBaffleEdgePoint() const
+{
+ return (type_ == vtInternalFeatureEdgeBaffle);
+}
template<class Gt, class Vb>
inline bool CGAL::indexedVertex<Gt, Vb>::externalBoundaryPoint() const
@@ -325,13 +326,15 @@ inline bool CGAL::indexedVertex<Gt, Vb>::externalBoundaryPoint() const
}
template<class Gt, class Vb>
-inline bool CGAL::indexedVertex<Gt, Vb>::externalBafflePoint() const
+inline bool CGAL::indexedVertex<Gt, Vb>::externalBaffleSurfacePoint() const
+{
+ return (type_ == vtExternalSurfaceBaffle);
+}
+
+template<class Gt, class Vb>
+inline bool CGAL::indexedVertex<Gt, Vb>::externalBaffleEdgePoint() const
{
- return
- (
- type_ == vtExternalSurfaceBaffle
- || type_ == vtExternalFeatureEdgeBaffle
- );
+ return (type_ == vtExternalFeatureEdgeBaffle);
}
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
index 43d9c770134a330f62629b67c488b32a91a1173d..ed944d86ce9642a35f9dcb5992a147c127947eef 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C
@@ -83,6 +83,13 @@ void Foam::conformationSurfaces::hasBoundedVolume
const pointField& surfPts = triSurf.points();
+ Info<< " Checking " << surface.name() << endl;
+
+ label nBaffles = 0;
+
+ Info<< " Index = " << surfaces_[s] << endl;
+ Info<< " Offset = " << regionOffset_[s] << endl;
+
forAll(triSurf, sI)
{
const label patchID =
@@ -98,7 +105,13 @@ void Foam::conformationSurfaces::hasBoundedVolume
{
sum += triSurf[sI].normal(surfPts);
}
+ else
+ {
+ nBaffles++;
+ }
}
+ Info<< " has " << nBaffles << " baffles out of "
+ << triSurf.size() << " triangles" << endl;
totalTriangles += triSurf.size();
}
@@ -153,8 +166,9 @@ void Foam::conformationSurfaces::readFeatures
{
const searchableSurface& surface = allGeometry_[surfaces_[surfI]];
- Info<< " features: " << surface.name() << " of type "
- << surface.type() << endl;
+ Info<< " features: " << surface.name()
+ << " of type " << surface.type()
+ << ", id: " << featureIndex << endl;
autoPtr<searchableSurfaceFeatures> ssFeatures
(
@@ -210,7 +224,8 @@ void Foam::conformationSurfaces::readFeatures
{
fileName feMeshName(featureDict.lookup("extendedFeatureEdgeMesh"));
- Info<< " features: " << feMeshName << endl;
+ Info<< " features: " << feMeshName << ", id: " << featureIndex
+ << endl;
features_.set
(
@@ -716,7 +731,11 @@ Foam::Field<bool> Foam::conformationSurfaces::wellInOutSide
const searchableSurface& surface(allGeometry_[surfaces_[s]]);
- if (!surface.hasVolumeType())
+ if
+ (
+ !surface.hasVolumeType()
+ //&& surfaceVolumeTests[s][i] == volumeType::UNKNOWN
+ )
{
pointField sample(1, samplePts[i]);
scalarField nearestDistSqr(1, GREAT);
@@ -733,7 +752,7 @@ Foam::Field<bool> Foam::conformationSurfaces::wellInOutSide
findSurfaceNearestIntersection
(
samplePts[i],
- info[0].rawPoint(),
+ info[0].rawPoint() - 1e-3*mag(hitDir)*hitDir,
surfHit,
hitSurface
);
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/autoDensity/autoDensity.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/autoDensity/autoDensity.C
index 5649e61480b80a056a7ec68e0054bd0453ee063f..3552323d68000af903fa0e671f6548f41217bedb 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/autoDensity/autoDensity.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/initialPointsMethod/autoDensity/autoDensity.C
@@ -187,7 +187,7 @@ Foam::label Foam::autoDensity::recurseAndFill
word recursionName
) const
{
- label maxDepth = 0;
+ label treeDepth = 0;
for (direction i = 0; i < 8; i++)
{
@@ -201,10 +201,10 @@ Foam::label Foam::autoDensity::recurseAndFill
{
if (levelLimit > 0)
{
- maxDepth =
+ treeDepth =
max
(
- maxDepth,
+ treeDepth,
recurseAndFill
(
initialPoints,
@@ -229,10 +229,10 @@ Foam::label Foam::autoDensity::recurseAndFill
if (!fillBox(initialPoints, subBB, true))
{
- maxDepth =
+ treeDepth =
max
(
- maxDepth,
+ treeDepth,
recurseAndFill
(
initialPoints,
@@ -259,10 +259,10 @@ Foam::label Foam::autoDensity::recurseAndFill
if (!fillBox(initialPoints, subBB, false))
{
- maxDepth =
+ treeDepth =
max
(
- maxDepth,
+ treeDepth,
recurseAndFill
(
initialPoints,
@@ -286,7 +286,7 @@ Foam::label Foam::autoDensity::recurseAndFill
}
}
- return maxDepth + 1;
+ return treeDepth + 1;
}
@@ -941,7 +941,7 @@ List<Vb::Point> autoDensity::initialPoints() const
Pout<< " Filling box " << hierBB << endl;
}
- label maxDepth = recurseAndFill
+ label treeDepth = recurseAndFill
(
initialPoints,
hierBB,
@@ -966,7 +966,7 @@ List<Vb::Point> autoDensity::initialPoints() const
<< scalar(nInitialPoints)/scalar(max(globalTrialPoints_, 1))
<< " success rate" << nl
<< indent
- << returnReduce(maxDepth, maxOp<label>())
+ << returnReduce(treeDepth, maxOp<label>())
<< " levels of recursion (maximum)"
<< decrIndent << decrIndent
<< endl;
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.C
new file mode 100644
index 0000000000000000000000000000000000000000..901d81ab3371d830ecda2098955e31831f27c4bd
--- /dev/null
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.C
@@ -0,0 +1,244 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "pointPairs.H"
+
+// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
+
+template<class Triangulation>
+inline Foam::Pair<Foam::labelPair>
+Foam::pointPairs<Triangulation>::orderPointPair
+(
+ const labelPair& vA,
+ const labelPair& vB
+) const
+{
+ return
+ (
+ (vA < vB)
+ ? Pair<labelPair>(vA, vB)
+ : Pair<labelPair>(vB, vA)
+ );
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::findPointPair
+(
+ const labelPair& vA,
+ const labelPair& vB
+) const
+{
+ if (vA == vB)
+ {
+ return false;
+ }
+ else if (find(orderPointPair(vA, vB)) == end())
+ {
+ return false;
+ }
+
+ return true;
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::insertPointPair
+(
+ const labelPair& vA,
+ const labelPair& vB
+)
+{
+ if (vA == vB)
+ {
+ return false;
+ }
+
+ return insert(orderPointPair(vA, vB));
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class Triangulation>
+Foam::pointPairs<Triangulation>::pointPairs(const Triangulation& triangulation)
+:
+ ptPairTable(),
+ triangulation_(triangulation)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class Triangulation>
+Foam::pointPairs<Triangulation>::~pointPairs()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::addPointPair
+(
+ const labelPair& vA,
+ const labelPair& vB
+)
+{
+ return insertPointPair(vA, vB);
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::addPointPair
+(
+ const labelPair& master,
+ const DynamicList<labelPair>& slaves
+)
+{
+ forAll(slaves, sI)
+ {
+ const labelPair& slave = slaves[sI];
+
+ addPointPair(master, slave);
+ }
+
+ return true;
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::addPointPair
+(
+ const label vA,
+ const label vB
+)
+{
+ const label procNo = Pstream::myProcNo();
+
+ labelPair a(vA, procNo);
+ labelPair b(vB, procNo);
+
+ return addPointPair(a, b);
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::isPointPair
+(
+ const Vertex_handle& vA,
+ const Vertex_handle& vB
+) const
+{
+ if (vA->boundaryPoint() && vB->boundaryPoint())
+ {
+ labelPair a(vA->index(), vA->procIndex());
+ labelPair b(vB->index(), vB->procIndex());
+
+ return findPointPair(a, b);
+ }
+
+ return false;
+}
+
+
+template<class Triangulation>
+inline bool Foam::pointPairs<Triangulation>::isPointPair
+(
+ const labelPair& vA,
+ const labelPair& vB
+) const
+{
+ return findPointPair(vA, vB);
+}
+
+
+template<class Triangulation>
+void Foam::pointPairs<Triangulation>::reIndex(const Map<label>& oldToNewIndices)
+{
+ pointPairs<Triangulation> newTable(triangulation_);
+
+ forAllConstIter(pointPairs, *this, iter)
+ {
+ Pair<labelPair> e = iter.key();
+
+ labelPair& start = e.first();
+ labelPair& end = e.second();
+
+ bool insert = true;
+
+ if (start.second() == Pstream::myProcNo())
+ {
+ Map<label>::const_iterator iter2 =
+ oldToNewIndices.find(start.first());
+
+ if (iter2 != oldToNewIndices.end())
+ {
+ if (iter2() != -1)
+ {
+ start.first() = iter2();
+ }
+ else
+ {
+ insert = false;
+ }
+ }
+ }
+
+ if (end.second() == Pstream::myProcNo())
+ {
+ Map<label>::const_iterator iter2 =
+ oldToNewIndices.find(end.first());
+
+ if (iter2 != oldToNewIndices.end())
+ {
+ if (iter2() != -1)
+ {
+ end.first() = iter2();
+ }
+ else
+ {
+ insert = false;
+ }
+ }
+ }
+
+ if (insert)
+ {
+ if (e.first() < e.second())
+ {
+ newTable.insert(e);
+ }
+ else if (e.first() > e.second())
+ {
+ newTable.insert(reverse(e));
+ }
+ }
+ }
+
+ this->transfer(newTable);
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.H
new file mode 100644
index 0000000000000000000000000000000000000000..b5fcd413eea867c3477e46c04261e0849e2f6766
--- /dev/null
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/pointPairs/pointPairs.H
@@ -0,0 +1,163 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ pointPairs
+
+Description
+ HashSet of unique edges. The edges are stored as a pair of pairs:
+
+ ( (local index, processor index) (local index, processor index) )
+
+ e.g.,
+
+ ( (0 1) (3 1) )
+ ( (0 2) (5 1) )
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef pointPairs_H
+#define pointPairs_H
+
+#include "labelPair.H"
+#include "HashSet.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+typedef HashSet
+<
+ Pair<labelPair>,
+ FixedList<labelPair, 2>::Hash<>
+> ptPairTable;
+
+/*---------------------------------------------------------------------------*\
+ Class pointPairs Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Triangulation>
+class pointPairs
+:
+ public ptPairTable
+{
+ // Private typedefs
+
+ typedef typename Triangulation::Vertex_handle Vertex_handle;
+
+
+ // Private data
+
+ const Triangulation& triangulation_;
+
+
+ // Private Member Functions
+
+ inline Pair<labelPair> orderPointPair
+ (
+ const labelPair& vA,
+ const labelPair& vB
+ ) const;
+
+ inline bool insertPointPair
+ (
+ const labelPair& vA,
+ const labelPair& vB
+ );
+
+ inline bool findPointPair
+ (
+ const labelPair& vA,
+ const labelPair& vB
+ ) const;
+
+
+public:
+
+ // Constructors
+
+ //- Construct from triangulation
+ pointPairs(const Triangulation& triangulation);
+
+
+ //- Destructor
+ ~pointPairs();
+
+
+ // Member Functions
+
+ // Access
+
+ inline bool isPointPair
+ (
+ const Vertex_handle& vA,
+ const Vertex_handle& vB
+ ) const;
+
+ inline bool isPointPair
+ (
+ const labelPair& vA,
+ const labelPair& vB
+ ) const;
+
+
+ // Edit
+
+ inline bool addPointPair
+ (
+ const labelPair& vA,
+ const labelPair& vB
+ );
+
+ inline bool addPointPair
+ (
+ const labelPair& master,
+ const DynamicList<labelPair>& slaves
+ );
+
+ inline bool addPointPair
+ (
+ const label vA,
+ const label vB
+ );
+
+ void reIndex(const Map<label>& oldToNewIndices);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+# include "pointPairs.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchableBoxFeatures.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchableBoxFeatures.C
index 78149fac1669328e7406e5cdfe888cd92f6d3e9a..76726ccc07350008afd8393e91d08b0ac3d8291d 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchableBoxFeatures.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchableBoxFeatures.C
@@ -104,11 +104,14 @@ Foam::searchableBoxFeatures::features() const
edgeNormals[10][0] = 1; edgeNormals[10][1] = 3;
edgeNormals[11][0] = 3; edgeNormals[11][1] = 0;
+ tmp<pointField> surfacePointsTmp(surface().points());
+ pointField& surfacePoints = surfacePointsTmp();
+
forAll(edgeDirections, eI)
{
edgeDirections[eI] =
- surface().points()()[treeBoundBox::edges[eI].end()]
- - surface().points()()[treeBoundBox::edges[eI].start()];
+ surfacePoints[treeBoundBox::edges[eI].end()]
+ - surfacePoints[treeBoundBox::edges[eI].start()];
normalDirections[eI] = labelList(2, 0);
for (label j = 0; j < 2; ++j)
@@ -116,8 +119,8 @@ Foam::searchableBoxFeatures::features() const
const vector cross =
(faceNormals[edgeNormals[eI][j]] ^ edgeDirections[eI]);
const vector fC0tofE0 =
- 0.5*(max(surface().points()() + min(surface().points()())))
- - surface().points()()[treeBoundBox::edges[eI].start()];
+ 0.5*(max(surfacePoints + min(surfacePoints)))
+ - surfacePoints[treeBoundBox::edges[eI].start()];
normalDirections[eI][j] =
(
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.C
new file mode 100644
index 0000000000000000000000000000000000000000..a61862c43211820cdff1a05d0494e429770fb325
--- /dev/null
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.C
@@ -0,0 +1,218 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "searchablePlateFeatures.H"
+#include "addToRunTimeSelectionTable.H"
+#include "treeBoundBox.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+defineTypeNameAndDebug(searchablePlateFeatures, 0);
+addToRunTimeSelectionTable
+(
+ searchableSurfaceFeatures,
+ searchablePlateFeatures,
+ dict
+);
+
+//! \cond - skip documentation : local scope only
+const Foam::label edgesArray[4][2] =
+{
+ {0, 1}, // 0
+ {0, 3},
+ {2, 1}, // 2
+ {2, 3}
+};
+//! \endcond
+
+const edgeList searchablePlateFeatures::edges(calcEdges(edgesArray));
+
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+Foam::edgeList Foam::searchablePlateFeatures::calcEdges
+(
+ const label edgesArray[4][2]
+)
+{
+ edgeList edges(4);
+ forAll(edges, edgeI)
+ {
+ edges[edgeI][0] = edgesArray[edgeI][0];
+ edges[edgeI][1] = edgesArray[edgeI][1];
+ }
+ return edges;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::searchablePlateFeatures::searchablePlateFeatures
+(
+ const searchableSurface& surface,
+ const dictionary& dict
+)
+:
+ searchableSurfaceFeatures(surface, dict),
+ mode_
+ (
+ extendedFeatureEdgeMesh::sideVolumeTypeNames_
+ [
+ dict.lookupOrDefault<word>("meshableSide", "inside")
+ ]
+ )
+{
+ Info<< indent
+ << " Meshable region = "
+ << extendedFeatureEdgeMesh::sideVolumeTypeNames_[mode_]
+ << endl;
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::searchablePlateFeatures::~searchablePlateFeatures()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::extendedFeatureEdgeMesh>
+Foam::searchablePlateFeatures::features() const
+{
+ autoPtr<extendedFeatureEdgeMesh> features;
+
+ vectorField faceNormals(1);
+ surface().getNormal(List<pointIndexHit>(1), faceNormals);
+
+ vectorField edgeDirections(4);
+ labelListList normalDirections(4);
+
+ labelListList edgeNormals(4);
+ forAll(edgeNormals, eI)
+ {
+ edgeNormals[eI].setSize(2, 0);
+ }
+ edgeNormals[0][0] = 0; edgeNormals[0][1] = 0;
+ edgeNormals[1][0] = 0; edgeNormals[1][1] = 0;
+ edgeNormals[2][0] = 0; edgeNormals[2][1] = 0;
+ edgeNormals[3][0] = 0; edgeNormals[3][1] = 0;
+
+ forAll(edgeDirections, eI)
+ {
+ edgeDirections[eI] =
+ surface().points()()[edges[eI].end()]
+ - surface().points()()[edges[eI].start()];
+
+ normalDirections[eI] = labelList(2, 0);
+ for (label j = 0; j < 2; ++j)
+ {
+ const vector cross =
+ (faceNormals[edgeNormals[eI][j]] ^ edgeDirections[eI]);
+ const vector fC0tofE0 =
+ 0.5*(max(surface().points()() + min(surface().points()())))
+ - surface().points()()[edges[eI].start()];
+
+ normalDirections[eI][j] =
+ (
+ (
+ (cross/(mag(cross) + VSMALL))
+ & (fC0tofE0/(mag(fC0tofE0)+ VSMALL))
+ )
+ > 0.0
+ ? 1
+ : -1
+ );
+ }
+ }
+
+ labelListList featurePointNormals(4);
+ labelListList featurePointEdges(4);
+ forAll(featurePointNormals, pI)
+ {
+ labelList& ftPtEdges = featurePointEdges[pI];
+ ftPtEdges.setSize(2, 0);
+
+ label edgeI = 0;
+ forAll(edges, eI)
+ {
+ const edge& e = edges[eI];
+
+ if (e.start() == pI)
+ {
+ ftPtEdges[edgeI++] = eI;
+ }
+ else if (e.end() == pI)
+ {
+ ftPtEdges[edgeI++] = eI;
+ }
+ }
+
+ labelList& ftPtNormals = featurePointNormals[pI];
+ ftPtNormals.setSize(1, 0);
+
+ ftPtNormals[0] = edgeNormals[ftPtEdges[0]][0];
+ }
+
+ labelList regionEdges;
+
+ features.set
+ (
+ new extendedFeatureEdgeMesh
+ (
+ IOobject
+ (
+ surface().name(),
+ surface().instance(),
+ "extendedFeatureEdgeMesh",
+ surface().db(),
+ IOobject::NO_READ,
+ IOobject::AUTO_WRITE
+ ),
+ surface().points(),
+ edges,
+ 4, 4, 4,
+ 0, 0, 4, 4, // 4 flat edges
+ faceNormals,
+ List<extendedFeatureEdgeMesh::sideVolumeType>(4, mode_),
+ edgeDirections,
+ normalDirections,
+ edgeNormals,
+ featurePointNormals,
+ featurePointEdges,
+ regionEdges
+ )
+ );
+
+ return features;
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.H
new file mode 100644
index 0000000000000000000000000000000000000000..3854a67d1fa427d305f7c7be1d0c4395ff70cdff
--- /dev/null
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/searchablePlateFeatures.H
@@ -0,0 +1,120 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::searchablePlateFeatures
+
+Description
+
+SourceFiles
+ searchablePlateFeatures.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef searchablePlateFeatures_H
+#define searchablePlateFeatures_H
+
+#include "searchableSurfaceFeatures.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class searchablePlateFeatures Declaration
+\*---------------------------------------------------------------------------*/
+
+class searchablePlateFeatures
+:
+ public searchableSurfaceFeatures
+{
+private:
+
+ //- To initialise edges.
+ static edgeList calcEdges(const label[4][2]);
+
+
+ // Private Member Data
+
+ //- Which side of the plate to mesh
+ const extendedFeatureEdgeMesh::sideVolumeType mode_;
+
+
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ searchablePlateFeatures(const searchablePlateFeatures&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const searchablePlateFeatures&);
+
+
+public:
+
+ //- Runtime type information
+ TypeName("searchablePlateFeatures");
+
+
+ // Static data members
+
+ //- Edge to point addressing
+ static const edgeList edges;
+
+
+ // Constructors
+
+ //- Construct from searchable surface and dictionary
+ searchablePlateFeatures
+ (
+ const searchableSurface& surface,
+ const dictionary& dict
+ );
+
+
+ //- Destructor
+ virtual ~searchablePlateFeatures();
+
+
+ // Member Functions
+
+ //- Return true for a searchable plate having features
+ virtual bool hasFeatures() const
+ {
+ return true;
+ }
+
+ //- Return an extendedFeatureEdgeMesh containing the features
+ virtual autoPtr<extendedFeatureEdgeMesh> features() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.C b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.C
index 3af1dd80d8dd80dc913618f284da2976ba6244c6..c61b92f6f1eb92481f936195df6cfcdfaf8ec265 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.C
@@ -57,10 +57,19 @@ Foam::triSurfaceMeshFeatures::triSurfaceMeshFeatures
)
:
searchableSurfaceFeatures(surface, dict),
- includedAngle_(readScalar(dict.lookup("includedAngle")))
+ includedAngle_(readScalar(dict.lookup("includedAngle"))),
+ mode_
+ (
+ extendedFeatureEdgeMesh::sideVolumeTypeNames_
+ [
+ dict.lookupOrDefault<word>("meshableSide", "inside")
+ ]
+ )
{
Info<< indent
- << " Included angle = " << includedAngle_
+ << " Included angle = " << includedAngle_ << nl
+ << " Meshable region = "
+ << extendedFeatureEdgeMesh::sideVolumeTypeNames_[mode_]
<< endl;
}
@@ -82,7 +91,12 @@ Foam::triSurfaceMeshFeatures::features() const
surfaceFeatures sFeat(surfMesh, includedAngle_);
- boolList surfBaffleRegions(surfMesh.patches().size(), false);
+ // @todo Need to read on a per region basis
+ boolList surfBaffleRegions
+ (
+ surfMesh.patches().size(),
+ (mode_ == extendedFeatureEdgeMesh::BOTH ? true : false)
+ );
features.set
(
@@ -90,8 +104,8 @@ Foam::triSurfaceMeshFeatures::features() const
(
sFeat,
surface().db(),
- surface().name() + ".extendedFeatureEdgeMesh"//,
- //surfBaffleRegions
+ surface().name() + ".extendedFeatureEdgeMesh",
+ surfBaffleRegions
)
);
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.H b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.H
index 190b7859ceac76d4d0c35cf402bafc50d2322e06..3105a38943b15f51b687b08edeeb626a6eaf6460 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.H
+++ b/applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/searchableSurfaceFeatures/triSurfaceMeshFeatures.H
@@ -55,6 +55,9 @@ private:
const scalar includedAngle_;
+ //- Which side of the surface to mesh
+ const extendedFeatureEdgeMesh::sideVolumeType mode_;
+
// Private Member Functions
diff --git a/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C b/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C
index f435aa616a7e4f8ef1c682d8dff461dc4e6494d5..c27d816cacfcee2b89b0ec32b1dd0d3c025ac446 100644
--- a/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C
+++ b/applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C
@@ -54,8 +54,6 @@ int main(int argc, char *argv[])
"conform to the initial points without any point motion"
);
- #include "addOverwriteOption.H"
-
#include "setRootCase.H"
#include "createTime.H"
@@ -63,7 +61,6 @@ int main(int argc, char *argv[])
const bool checkGeometry = args.optionFound("checkGeometry");
const bool conformationOnly = args.optionFound("conformationOnly");
- const bool overwrite = args.optionFound("overwrite");
IOdictionary foamyHexMeshDict
(
@@ -123,10 +120,7 @@ int main(int argc, char *argv[])
{
mesh.initialiseForConformation();
- if (!overwrite)
- {
- runTime++;
- }
+ runTime++;
mesh.writeMesh(runTime.timeName());
}
@@ -134,8 +128,10 @@ int main(int argc, char *argv[])
{
mesh.initialiseForMotion();
- while (runTime.loop())
+ while (runTime.run())
{
+ runTime++;
+
Info<< nl << "Time = " << runTime.timeName() << endl;
mesh.move();
diff --git a/applications/utilities/mesh/generation/foamyQuadMesh/Make/options b/applications/utilities/mesh/generation/foamyQuadMesh/Make/options
index 62aec78cdb2e485624764ac37477d59dd9387207..7c52fa131a3661b0e93237c8e67c655406a7f8a4 100755
--- a/applications/utilities/mesh/generation/foamyQuadMesh/Make/options
+++ b/applications/utilities/mesh/generation/foamyQuadMesh/Make/options
@@ -9,6 +9,7 @@ EXE_INC = \
${EXE_FROUNDING_MATH} \
${EXE_NDEBUG} \
${CGAL_INC} \
+ ${c++CGALWARN} \
-I$(FOAM_APP)/utilities/mesh/generation/extrude2DMesh/extrude2DMesh/lnInclude \
-I../cvMesh/vectorTools \
-IconformalVoronoi2DMesh/lnInclude \
diff --git a/applications/utilities/mesh/manipulation/orientFaceZone/Make/files b/applications/utilities/mesh/manipulation/orientFaceZone/Make/files
index 07214af76c6b92c6f0ba9050670b24a5c35cb248..f5a0821f85e6c216410c83854236e9975a46c698 100644
--- a/applications/utilities/mesh/manipulation/orientFaceZone/Make/files
+++ b/applications/utilities/mesh/manipulation/orientFaceZone/Make/files
@@ -1,4 +1,3 @@
-patchFaceOrientation.C
orientFaceZone.C
EXE = $(FOAM_APPBIN)/orientFaceZone
diff --git a/applications/utilities/mesh/manipulation/orientFaceZone/Make/options b/applications/utilities/mesh/manipulation/orientFaceZone/Make/options
index 9f08e8d2a80f6f1fb269c1815170453ae49cf783..30fde688e5aa0fe88f1940882bf0b5d3f888af9c 100644
--- a/applications/utilities/mesh/manipulation/orientFaceZone/Make/options
+++ b/applications/utilities/mesh/manipulation/orientFaceZone/Make/options
@@ -1,7 +1,9 @@
EXE_INC = \
-I$(LIB_SRC)/meshTools/lnInclude \
+ -I$(LIB_SRC)/mesh/autoMesh/lnInclude \
-I$(LIB_SRC)/triSurface/lnInclude
EXE_LIBS = \
-lmeshTools \
+ -lautoMesh \
-ltriSurface
diff --git a/applications/utilities/surface/surfaceBooleanFeatures/surfaceBooleanFeatures.C b/applications/utilities/surface/surfaceBooleanFeatures/surfaceBooleanFeatures.C
index debb7b05887894fd43e7a6a2f2843873705df512..e5c928e29889659ac60ad3220e20a4d261a524c0 100644
--- a/applications/utilities/surface/surfaceBooleanFeatures/surfaceBooleanFeatures.C
+++ b/applications/utilities/surface/surfaceBooleanFeatures/surfaceBooleanFeatures.C
@@ -78,11 +78,79 @@ Description
#include "booleanSurface.H"
#include "edgeIntersections.H"
#include "meshTools.H"
+#include "labelPair.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+bool intersectSurfaces
+(
+ triSurface& surf,
+ edgeIntersections& edgeCuts
+)
+{
+ bool hasMoved = false;
+
+ for (label iter = 0; iter < 10; iter++)
+ {
+ Info<< "Determining intersections of surface edges with itself" << endl;
+
+ // Determine surface edge intersections. Allow surface to be moved.
+
+ // Number of iterations needed to resolve degenerates
+ label nIters = 0;
+ {
+ triSurfaceSearch querySurf(surf);
+
+ scalarField surfPointTol
+ (
+ 1e-3*edgeIntersections::minEdgeLength(surf)
+ );
+
+ // Determine raw intersections
+ edgeCuts = edgeIntersections
+ (
+ surf,
+ querySurf,
+ surfPointTol
+ );
+
+ // Shuffle a bit to resolve degenerate edge-face hits
+ {
+ pointField points(surf.points());
+
+ nIters =
+ edgeCuts.removeDegenerates
+ (
+ 5, // max iterations
+ surf,
+ querySurf,
+ surfPointTol,
+ points // work array
+ );
+
+ if (nIters != 0)
+ {
+ // Update geometric quantities
+ surf.movePoints(points);
+ hasMoved = true;
+ }
+ }
+ }
+ }
+
+ if (hasMoved)
+ {
+ fileName newFile("surf.obj");
+ Info<< "Surface has been moved. Writing to " << newFile << endl;
+ surf.write(newFile);
+ }
+
+ return hasMoved;
+}
+
+
// Keep on shuffling surface points until no more degenerate intersections.
// Moves both surfaces and updates set of edge cuts.
bool intersectSurfaces
@@ -238,6 +306,106 @@ label calcNormalDirection
}
+void calcEdgeCuts
+(
+ triSurface& surf1,
+ triSurface& surf2,
+ const bool perturb,
+ edgeIntersections& edge1Cuts,
+ edgeIntersections& edge2Cuts
+)
+{
+ if (perturb)
+ {
+ intersectSurfaces
+ (
+ surf1,
+ edge1Cuts,
+ surf2,
+ edge2Cuts
+ );
+ }
+ else
+ {
+ triSurfaceSearch querySurf2(surf2);
+
+ Info<< "Determining intersections of surf1 edges with surf2 faces"
+ << endl;
+
+ edge1Cuts = edgeIntersections
+ (
+ surf1,
+ querySurf2,
+ 1e-3*edgeIntersections::minEdgeLength(surf1)
+ );
+
+ triSurfaceSearch querySurf1(surf1);
+
+ Info<< "Determining intersections of surf2 edges with surf1 faces"
+ << endl;
+
+ edge2Cuts = edgeIntersections
+ (
+ surf2,
+ querySurf1,
+ 1e-3*edgeIntersections::minEdgeLength(surf2)
+ );
+ }
+}
+
+
+void calcFeaturePoints(const pointField& points, const edgeList& edges)
+{
+ edgeMesh eMesh(points, edges);
+
+ const labelListList& pointEdges = eMesh.pointEdges();
+
+
+ // Get total number of feature points
+ label nFeaturePoints = 0;
+ forAll(pointEdges, pI)
+ {
+ const labelList& pEdges = pointEdges[pI];
+
+ if (pEdges.size() == 1)
+ {
+ nFeaturePoints++;
+ }
+ }
+
+
+ // Calculate addressing from feature point to cut point and cut edge
+ labelList featurePointToCutPoint(nFeaturePoints);
+ labelList featurePointToCutEdge(nFeaturePoints);
+
+ label nFeatPts = 0;
+ forAll(pointEdges, pI)
+ {
+ const labelList& pEdges = pointEdges[pI];
+
+ if (pEdges.size() == 1)
+ {
+ featurePointToCutPoint[nFeatPts] = pI;
+ featurePointToCutEdge[nFeatPts] = pEdges[0];
+ nFeatPts++;
+ }
+ }
+
+
+
+ label concaveStart = 0;
+ label mixedStart = 0;
+ label nonFeatureStart = nFeaturePoints;
+
+
+ labelListList featurePointNormals(nFeaturePoints);
+ labelListList featurePointEdges(nFeaturePoints);
+ labelList regionEdges;
+
+
+}
+
+
int main(int argc, char *argv[])
{
argList::noParallel();
@@ -319,52 +487,33 @@ int main(int argc, char *argv[])
<< exit(FatalError);
}
- if (args.optionFound("perturb"))
- {
- intersectSurfaces
- (
- surf1,
- edge1Cuts,
- surf2,
- edge2Cuts
- );
- }
- else
- {
- triSurfaceSearch querySurf2(surf2);
-
- Info<< "Determining intersections of surf1 edges with surf2 faces"
- << endl;
-
- edge1Cuts = edgeIntersections
- (
- surf1,
- querySurf2,
- 1e-3*edgeIntersections::minEdgeLength(surf1)
- );
-
- triSurfaceSearch querySurf1(surf1);
-
- Info<< "Determining intersections of surf2 edges with surf1 faces"
- << endl;
-
- edge2Cuts = edgeIntersections
- (
- surf2,
- querySurf1,
- 1e-3*edgeIntersections::minEdgeLength(surf2)
- );
- }
+ // Calculate the points where the edges are cut by the other surface
+ calcEdgeCuts
+ (
+ surf1,
+ surf2,
+ args.optionFound("perturb"),
+ edge1Cuts,
+ edge2Cuts
+ );
- // Determine intersection edges
- surfaceIntersection inter(surf1, edge1Cuts, surf2, edge2Cuts);
+ // Determine intersection edges from the edge cuts
+ surfaceIntersection inter
+ (
+ surf1,
+ edge1Cuts,
+ surf2,
+ edge2Cuts
+ );
- fileName sFeatFileName =
+ fileName sFeatFileName
+ (
surf1Name.lessExt().name()
+ "_"
+ surf2Name.lessExt().name()
+ "_"
- + action;
+ + action
+ );
label nFeatEds = inter.cutEdges().size();
@@ -393,6 +542,7 @@ int main(int argc, char *argv[])
edgeDirections[cutEdgeI] = fE.vec(inter.cutPoints());
normals.append(norm1);
+ eNormals.append(normals.size() - 1);
if (surf1Baffle)
{
@@ -416,9 +566,8 @@ int main(int argc, char *argv[])
);
}
- eNormals.append(normals.size() - 1);
-
normals.append(norm2);
+ eNormals.append(normals.size() - 1);
if (surf2Baffle)
{
@@ -443,8 +592,6 @@ int main(int argc, char *argv[])
);
}
- eNormals.append(normals.size() - 1);
-
if (surf1Baffle)
{
@@ -509,6 +656,38 @@ int main(int argc, char *argv[])
label internalStart = -1;
+ label nIntOrExt = 0;
+ label nFlat = 0;
+ label nOpen = 0;
+ label nMultiple = 0;
+
+ forAll(edgeNormals, eI)
+ {
+ label nEdNorms = edgeNormals[eI].size();
+
+ if (nEdNorms == 1)
+ {
+ nOpen++;
+ }
+ else if (nEdNorms == 2)
+ {
+ const vector& n0(normals[edgeNormals[eI][0]]);
+ const vector& n1(normals[edgeNormals[eI][1]]);
+
+ if ((n0 & n1) > extendedFeatureEdgeMesh::cosNormalAngleTol_)
+ {
+ nFlat++;
+ }
+ else
+ {
+ nIntOrExt++;
+ }
+ }
+ else if (nEdNorms > 2)
+ {
+ nMultiple++;
+ }
+ }
if (validActions[action] == booleanSurface::UNION)
{
@@ -520,7 +699,7 @@ int main(int argc, char *argv[])
else
{
// All edges are external
- internalStart = nFeatEds;
+ internalStart = nIntOrExt;
}
}
else if (validActions[action] == booleanSurface::INTERSECTION)
@@ -528,7 +707,7 @@ int main(int argc, char *argv[])
if (!invertedSpace)
{
// All edges are external
- internalStart = nFeatEds;
+ internalStart = nIntOrExt;
}
else
{
@@ -539,7 +718,7 @@ int main(int argc, char *argv[])
else if (validActions[action] == booleanSurface::DIFFERENCE)
{
// All edges are external
- internalStart = nFeatEds;
+ internalStart = nIntOrExt;
}
else
{
@@ -569,6 +748,8 @@ int main(int argc, char *argv[])
normalDirectionsTmp[i] = normalDirections[i];
}
+ calcFeaturePoints(inter.cutPoints(), inter.cutEdges());
+
extendedFeatureEdgeMesh feMesh
(
IOobject
@@ -582,18 +763,22 @@ int main(int argc, char *argv[])
),
inter.cutPoints(),
inter.cutEdges(),
+
0, // concaveStart,
0, // mixedStart,
0, // nonFeatureStart,
+
internalStart, // internalStart,
- nFeatEds, // flatStart,
- nFeatEds, // openStart,
- nFeatEds, // multipleStart,
+ nIntOrExt, // flatStart,
+ nIntOrExt + nFlat, // openStart,
+ nIntOrExt + nFlat + nOpen, // multipleStart,
+
normalsTmp,
normalVolumeTypesTmp,
edgeDirections,
normalDirectionsTmp,
edgeNormalsTmp,
+
labelListList(0), // featurePointNormals,
labelListList(0), // featurePointEdges,
labelList(0) // regionEdges
diff --git a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
index 58f321297b764acb9e06b86a8119e919c7bc9d68..77582934808a52c251a6bc18968efb545830175c 100644
--- a/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
+++ b/applications/utilities/surface/surfaceFeatureExtract/surfaceFeatureExtract.C
@@ -925,7 +925,7 @@ void writeStats(const extendedFeatureEdgeMesh& fem, Ostream& os)
<< " open edges : "
<< (fem.multipleStart()- fem.openStart()) << nl
<< " multiply connected : "
- << (fem.edges().size()- fem.multipleStart()) << nl;
+ << (fem.edges().size()- fem.multipleStart()) << endl;
}
@@ -1223,17 +1223,34 @@ int main(int argc, char *argv[])
<< " internal edges : " << newSet.nInternalEdges() << nl
<< endl;
- //if (writeObj)
- //{
- // newSet.writeObj("final");
- //}
+ if (writeObj)
+ {
+ newSet.writeObj("final");
+ }
+
+ boolList surfBaffleRegions(surf.patches().size(), false);
+
+ wordList surfBaffleNames;
+ surfaceDict.readIfPresent("baffles", surfBaffleNames);
+
+ forAll(surf.patches(), pI)
+ {
+ const word& name = surf.patches()[pI].name();
+
+ if (findIndex(surfBaffleNames, name) != -1)
+ {
+ Info<< "Adding baffle region " << name << endl;
+ surfBaffleRegions[pI] = true;
+ }
+ }
// Extracting and writing a extendedFeatureEdgeMesh
extendedFeatureEdgeMesh feMesh
(
newSet,
runTime,
- sFeatFileName + ".extendedFeatureEdgeMesh"
+ sFeatFileName + ".extendedFeatureEdgeMesh",
+ surfBaffleRegions
);
diff --git a/src/OpenFOAM/db/error/error.H b/src/OpenFOAM/db/error/error.H
index 8216f5f71fa0e223fc711e1a367ecba292951095..97a08b175c612c4f214b1325b6e62296bd4dc014 100644
--- a/src/OpenFOAM/db/error/error.H
+++ b/src/OpenFOAM/db/error/error.H
@@ -28,11 +28,11 @@ Description
Class to handle errors and exceptions in a simple, consistent stream-based
manner.
- The error class is globaly instantiated with a title string. Errors,
+ The error class is globally instantiated with a title string. Errors,
messages and other data are piped to the messageStream class in the
standard manner. Manipulators are supplied for exit and abort which may
- terminate the program or throw an exception depending of if the exception
- handling has beed switched on (off by default).
+ terminate the program or throw an exception depending on whether the
+ exception handling has been switched on (off by default).
Usage
\code
diff --git a/src/OpenFOAM/memory/autoPtr/autoPtrI.H b/src/OpenFOAM/memory/autoPtr/autoPtrI.H
index 55e032fcd4548cebe3e8d99164ffc37e22d0e397..3557ffc49b265deadc351e5c83d931d1e19069d5 100644
--- a/src/OpenFOAM/memory/autoPtr/autoPtrI.H
+++ b/src/OpenFOAM/memory/autoPtr/autoPtrI.H
@@ -52,10 +52,14 @@ inline Foam::autoPtr<T>::autoPtr(const autoPtr<T>& ap, const bool reUse)
ptr_ = ap.ptr_;
ap.ptr_ = 0;
}
- else
+ else if (ap.valid())
{
ptr_ = ap().clone().ptr();
}
+ else
+ {
+ ptr_ = NULL;
+ }
}
diff --git a/src/dynamicMesh/polyMeshFilter/polyMeshFilter.C b/src/dynamicMesh/polyMeshFilter/polyMeshFilter.C
index 2958356dc96befdded8b23d9de6de809b5a1939c..968ea6c483595cab695831bdb4ba12e3a687a139 100644
--- a/src/dynamicMesh/polyMeshFilter/polyMeshFilter.C
+++ b/src/dynamicMesh/polyMeshFilter/polyMeshFilter.C
@@ -32,7 +32,9 @@ License
#include "polyTopoChange.H"
#include "globalIndex.H"
#include "PackedBoolList.H"
-#include "faceZone.H"
+#include "pointSet.H"
+#include "faceSet.H"
+#include "cellSet.H"
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
@@ -42,6 +44,26 @@ defineTypeNameAndDebug(polyMeshFilter, 0);
}
+void Foam::polyMeshFilter::updateSets(const mapPolyMesh& map)
+{
+ updateSets<pointSet>(map);
+ updateSets<faceSet>(map);
+ updateSets<cellSet>(map);
+}
+
+
+void Foam::polyMeshFilter::copySets
+(
+ const polyMesh& oldMesh,
+ const polyMesh& newMesh
+)
+{
+ copySets<pointSet>(oldMesh, newMesh);
+ copySets<faceSet>(oldMesh, newMesh);
+ copySets<cellSet>(oldMesh, newMesh);
+}
+
+
Foam::autoPtr<Foam::fvMesh> Foam::polyMeshFilter::copyMesh(const fvMesh& mesh)
{
polyTopoChange originalMeshToNewMesh(mesh);
@@ -72,6 +94,8 @@ Foam::autoPtr<Foam::fvMesh> Foam::polyMeshFilter::copyMesh(const fvMesh& mesh)
meshCopy().movePoints(map.preMotionPoints());
}
+ copySets(mesh, meshCopy());
+
return meshCopy;
}
@@ -359,6 +383,7 @@ Foam::label Foam::polyMeshFilter::filterFaces
{
newMesh.movePoints(newMap.preMotionPoints());
}
+ updateSets(newMap);
updatePointPriorities(newMesh, newMap.pointMap());
@@ -475,6 +500,7 @@ Foam::label Foam::polyMeshFilter::filterEdges
{
newMesh.movePoints(newMap.preMotionPoints());
}
+ updateSets(newMap);
// Synchronise the factors
mapOldMeshEdgeFieldToNewMesh
@@ -926,14 +952,14 @@ Foam::label Foam::polyMeshFilter::filter(const label nOriginalBadFaces)
}
-Foam::label Foam::polyMeshFilter::filter(const faceZone& fZone)
+Foam::label Foam::polyMeshFilter::filter(const faceSet& fSet)
{
minEdgeLen_.resize(mesh_.nEdges(), minLen());
faceFilterFactor_.resize(mesh_.nFaces(), initialFaceLengthFactor());
forAll(faceFilterFactor_, fI)
{
- if (fZone.whichFace(fI) == -1)
+ if (!fSet.found(fI))
{
faceFilterFactor_[fI] = -1;
}
diff --git a/src/dynamicMesh/polyMeshFilter/polyMeshFilter.H b/src/dynamicMesh/polyMeshFilter/polyMeshFilter.H
index 95cddd0058cf214b7ad9a46b8359b09bf128d6cf..816e75e126e8d1a66eb0dacab34f3fea6c862ba9 100644
--- a/src/dynamicMesh/polyMeshFilter/polyMeshFilter.H
+++ b/src/dynamicMesh/polyMeshFilter/polyMeshFilter.H
@@ -25,13 +25,14 @@ Class
Foam::polyMeshFilter
Description
- Filter the edges and faces of a polyMesh whilst satisfying the given mesh
+ Remove the edges and faces of a polyMesh whilst satisfying the given mesh
quality criteria.
Works on a copy of the mesh.
SourceFiles
polyMeshFilter.C
+ polyMeshFilterTemplates.C
\*---------------------------------------------------------------------------*/
@@ -53,7 +54,7 @@ namespace Foam
class polyMesh;
class fvMesh;
class PackedBoolList;
-class faceZone;
+class faceSet;
/*---------------------------------------------------------------------------*\
Class polyMeshFilter Declaration
@@ -71,8 +72,12 @@ class polyMeshFilter
//- Copy of the original mesh to perform the filtering on
autoPtr<fvMesh> newMeshPtr_;
- //-
+ //- Original point priorities. If a point has a higher priority than
+ // another point then the edge between them collapses towards the
+ // point with the higher priority. (e.g. 2 is a higher priority than 1)
labelList originalPointPriority_;
+
+ //- Point priority associated with the new mesh
autoPtr<labelList> pointPriority_;
//- The minimum edge length for each edge
@@ -84,6 +89,12 @@ class polyMeshFilter
// Private Member Functions
+ template<typename T>
+ static void updateSets(const mapPolyMesh& map);
+
+ template<typename T>
+ static void copySets(const polyMesh& oldMesh, const polyMesh& newMesh);
+
label filterFacesLoop(const label nOriginalBadFaces);
label filterFaces
@@ -209,6 +220,7 @@ public:
// mesh has actually been filtered.
const autoPtr<fvMesh>& filteredMesh() const;
+ //- Return the new pointPriority list.
const autoPtr<labelList>& pointPriority() const;
@@ -217,11 +229,21 @@ public:
//- Return a copy of an fvMesh
static autoPtr<fvMesh> copyMesh(const fvMesh& mesh);
+ //- Copy loaded topoSets from the old mesh to the new mesh
+ static void copySets
+ (
+ const polyMesh& oldMesh,
+ const polyMesh& newMesh
+ );
+
+ //- Update the loaded topoSets
+ static void updateSets(const mapPolyMesh& map);
+
//- 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 all faces that are in the face set
+ label filter(const faceSet& fSet);
//- Filter edges only.
label filterEdges(const label nOriginalBadFaces);
@@ -234,6 +256,12 @@ public:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+#ifdef NoRepository
+# include "polyMeshFilterTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
#endif
// ************************************************************************* //
diff --git a/src/dynamicMesh/polyMeshFilter/polyMeshFilterSettings.H b/src/dynamicMesh/polyMeshFilter/polyMeshFilterSettings.H
index e248ce80ad972be723f4847eb0766be45a51cede..93bc4bd3b911e00234ca9b83833a9cf55e21b4e9 100644
--- a/src/dynamicMesh/polyMeshFilter/polyMeshFilterSettings.H
+++ b/src/dynamicMesh/polyMeshFilter/polyMeshFilterSettings.H
@@ -29,6 +29,7 @@ Description
SourceFiles
polyMeshFilterSettings.C
+ polyMeshFilterSettingsI.H
\*---------------------------------------------------------------------------*/
diff --git a/src/dynamicMesh/polyMeshFilter/polyMeshFilterTemplates.C b/src/dynamicMesh/polyMeshFilter/polyMeshFilterTemplates.C
new file mode 100644
index 0000000000000000000000000000000000000000..08fa9d6efe8d9845da689754878b68d4350c28d8
--- /dev/null
+++ b/src/dynamicMesh/polyMeshFilter/polyMeshFilterTemplates.C
@@ -0,0 +1,105 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+License
+ This file is part of OpenFOAM.
+
+ OpenFOAM is free software: you can redistribute it and/or modify it
+ under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+ ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+ FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
+ for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "polyMeshFilter.H"
+#include "polyMesh.H"
+#include "mapPolyMesh.H"
+#include "IOobjectList.H"
+
+// * * * * * * * * * * * * * Public Member Functions * * * * * * * * * * * * //
+
+template<typename SetType>
+void Foam::polyMeshFilter::updateSets(const mapPolyMesh& map)
+{
+ HashTable<const SetType*> sets =
+ map.mesh().objectRegistry::lookupClass<const SetType>();
+
+ forAllIter(typename HashTable<const SetType*>, sets, iter)
+ {
+ SetType& set = const_cast<SetType&>(*iter());
+ set.updateMesh(map);
+ set.sync(map.mesh());
+ }
+
+ IOobjectList Objects
+ (
+ map.mesh().time(),
+ map.mesh().facesInstance(),
+ "polyMesh/sets"
+ );
+
+ IOobjectList fileSets(Objects.lookupClass(SetType::typeName));
+
+ forAllConstIter(IOobjectList, fileSets, iter)
+ {
+ if (!sets.found(iter.key()))
+ {
+ // Not in memory. Load it.
+ SetType set(*iter());
+ set.updateMesh(map);
+
+ set.write();
+ }
+ }
+}
+
+
+template<typename SetType>
+void Foam::polyMeshFilter::copySets
+(
+ const polyMesh& oldMesh,
+ const polyMesh& newMesh
+)
+{
+ HashTable<const SetType*> sets =
+ oldMesh.objectRegistry::lookupClass<const SetType>();
+
+ forAllConstIter(typename HashTable<const SetType*>, sets, iter)
+ {
+ const SetType& set = *iter();
+
+ if (newMesh.objectRegistry::foundObject<SetType>(set.name()))
+ {
+ const SetType& origSet =
+ newMesh.objectRegistry::lookupObject<SetType>(set.name());
+
+ const_cast<SetType&>(origSet) = set;
+ const_cast<SetType&>(origSet).sync(newMesh);
+ }
+ else
+ {
+ SetType* newSet
+ (
+ new SetType(newMesh, set.name(), set, set.writeOpt())
+ );
+
+ newSet->store();
+ newSet->sync(newMesh);
+ }
+ }
+}
+
+
+// ************************************************************************* //
diff --git a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.C b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.C
index 54ad050fc65d61a433e58d252cad002261320c48..0d03912ceed215ebfc9fb6b96d267912bac487be 100644
--- a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.C
+++ b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.C
@@ -263,7 +263,8 @@ Foam::extendedFeatureEdgeMesh::extendedFeatureEdgeMesh
(
const surfaceFeatures& sFeat,
const objectRegistry& obr,
- const fileName& sFeatFileName
+ const fileName& sFeatFileName,
+ const boolList& surfBaffleRegions
)
:
regIOobject
@@ -337,7 +338,12 @@ Foam::extendedFeatureEdgeMesh::extendedFeatureEdgeMesh
// Check to see if the points have been already used
if (faceMap[eFI] == -1)
{
- normalVolumeTypes_[nAdded++] = INSIDE;
+ normalVolumeTypes_[nAdded++] =
+ (
+ surfBaffleRegions[surf[eFI].region()]
+ ? BOTH
+ : INSIDE
+ );
faceMap[eFI] = nAdded - 1;
}
@@ -492,7 +498,7 @@ Foam::extendedFeatureEdgeMesh::classifyEdge
const List<vector>& norms,
const labelList& edNorms,
const vector& fC0tofC1
-) const
+)
{
label nEdNorms = edNorms.size();
@@ -502,8 +508,8 @@ Foam::extendedFeatureEdgeMesh::classifyEdge
}
else if (nEdNorms == 2)
{
- const vector n0(norms[edNorms[0]]);
- const vector n1(norms[edNorms[1]]);
+ const vector& n0(norms[edNorms[0]]);
+ const vector& n1(norms[edNorms[1]]);
if ((n0 & n1) > cosNormalAngleTol_)
{
@@ -1342,6 +1348,22 @@ void Foam::extendedFeatureEdgeMesh::writeObj
meshTools::writeOBJ(regionStr, points()[e[1]]); verti++;
regionStr << "l " << verti-1 << ' ' << verti << endl;
}
+
+ OFstream edgeDirsStr(prefix + "_edgeDirections.obj");
+ Info<< "Writing edge directions to " << edgeDirsStr.name() << endl;
+
+ forAll(edgeDirections_, i)
+ {
+ const vector& eVec = edgeDirections_[i];
+ const edge& e = edges()[i];
+
+ meshTools::writeOBJ
+ (
+ edgeDirsStr,
+ points()[e.start()],
+ eVec + points()[e.start()]
+ );
+ }
}
diff --git a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.H b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.H
index dec5c416965d8c2a98bed81581f3f1a8ed67737f..3695a3901261e44256bb4c9c90474b94283a7c2d 100644
--- a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.H
+++ b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMesh.H
@@ -119,14 +119,14 @@ public:
static const Foam::NamedEnum<sideVolumeType, 4> sideVolumeTypeNames_;
+ //- Angular closeness tolerance for treating normals as the same
+ static scalar cosNormalAngleTol_;
+
private:
// Static data
- //- Angular closeness tolerance for treating normals as the same
- static scalar cosNormalAngleTol_;
-
//- Index of the start of the convex feature points - static as 0
static label convexStart_;
@@ -202,15 +202,6 @@ private:
// data for edges and normals.
pointStatus classifyFeaturePoint(label ptI) const;
- //- Classify the type of feature edge. Requires face centre 0 to face
- // centre 1 vector to distinguish internal from external
- edgeStatus classifyEdge
- (
- const List<vector>& norms,
- const labelList& edNorms,
- const vector& fC0tofC1
- ) const;
-
template<class Patch>
void sortPointsAndEdges
(
@@ -259,7 +250,8 @@ public:
(
const surfaceFeatures& sFeat,
const objectRegistry& obr,
- const fileName& sFeatFileName
+ const fileName& sFeatFileName,
+ const boolList& surfBaffleRegions
);
//- Construct from PrimitivePatch
@@ -434,6 +426,9 @@ public:
//- Return the edgeStatus of a specified edge
inline edgeStatus getEdgeStatus(label edgeI) const;
+ //- Return the baffle faces of a specified edge
+ inline PackedList<2> edgeBaffles(label edgeI) const;
+
//- Demand driven construction of octree for feature points
const indexedOctree<treeDataPoint>& pointTree() const;
@@ -468,6 +463,16 @@ public:
friend Istream& operator>>(Istream& is, sideVolumeType& vt);
friend Ostream& operator<<(Ostream& os, const sideVolumeType& vt);
+
+
+ //- Classify the type of feature edge. Requires face centre 0 to face
+ // centre 1 vector to distinguish internal from external
+ static edgeStatus classifyEdge
+ (
+ const List<vector>& norms,
+ const labelList& edNorms,
+ const vector& fC0tofC1
+ );
};
diff --git a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshI.H b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshI.H
index a0ba24daaef63c0f4df630f93366db24e74cf10a..9259b06cb52c61cfb20625371f79107364fb9f6f 100644
--- a/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshI.H
+++ b/src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshI.H
@@ -264,4 +264,27 @@ Foam::extendedFeatureEdgeMesh::getEdgeStatus(label edgeI) const
}
+inline Foam::PackedList<2> Foam::extendedFeatureEdgeMesh::edgeBaffles
+(
+ label edgeI
+) const
+{
+ const labelList& eNormals = edgeNormals_[edgeI];
+
+ DynamicList<label> edgeBaffles(eNormals.size());
+
+ forAll(eNormals, enI)
+ {
+ const label normI = eNormals[enI];
+
+ if (normalVolumeTypes_[normI])
+ {
+ edgeBaffles.append(normI);
+ }
+ }
+
+ return PackedList<2>(edgeBaffles);
+}
+
+
// ************************************************************************* //
diff --git a/src/fvAgglomerationMethods/Allwmake b/src/fvAgglomerationMethods/Allwmake
index abbcf4664b0210cd9ac9216ca3ec2bc776960d74..4fc3ed6029541a50fe93dd1e917281a254f80ed8 100755
--- a/src/fvAgglomerationMethods/Allwmake
+++ b/src/fvAgglomerationMethods/Allwmake
@@ -12,4 +12,22 @@ fi
wmake $makeType pairPatchAgglomeration
+
+## get SCOTCH_VERSION, SCOTCH_ARCH_PATH
+#if settings=`$WM_PROJECT_DIR/bin/foamEtcFile config/scotch.sh`
+#then
+# . $settings
+# echo "using SCOTCH_ARCH_PATH=$SCOTCH_ARCH_PATH"
+#else
+# echo
+# echo "Error: no config/scotch.sh settings"
+# echo
+#fi
+#
+#if [ -n "$SCOTCH_ARCH_PATH" ]
+#then
+# wmake $makeType scotchGamgAgglomeration
+#fi
+
+
# ----------------------------------------------------------------- end-of-file
diff --git a/src/mesh/autoMesh/Make/files b/src/mesh/autoMesh/Make/files
index d5bfe6f751fd427461185d9d9e7229c3b6c13ae0..dd0bfd0391fd54ceaa11fc0ceeb185f3f572f625 100644
--- a/src/mesh/autoMesh/Make/files
+++ b/src/mesh/autoMesh/Make/files
@@ -17,6 +17,8 @@ $(autoHexMesh)/meshRefinement/meshRefinement.C
$(autoHexMesh)/meshRefinement/meshRefinementMerge.C
$(autoHexMesh)/meshRefinement/meshRefinementProblemCells.C
$(autoHexMesh)/meshRefinement/meshRefinementRefine.C
+$(autoHexMesh)/meshRefinement/patchFaceOrientation.C
+
$(autoHexMesh)/refinementFeatures/refinementFeatures.C
$(autoHexMesh)/refinementSurfaces/surfaceZonesInfo.C
$(autoHexMesh)/refinementSurfaces/refinementSurfaces.C
diff --git a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/refinementParameters/refinementParameters.C b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/refinementParameters/refinementParameters.C
index 0e2fcd07ce7b088dba1a648903f38474869d293c..f9676fcbe6c592f0d61c1d4d1c10a0f16acb5515 100644
--- a/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/refinementParameters/refinementParameters.C
+++ b/src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/refinementParameters/refinementParameters.C
@@ -112,8 +112,6 @@ const
if (localCellI != -1)
{
- Pout<< "Found point " << keepPoint << " in cell " << localCellI
- << " on processor " << Pstream::myProcNo() << endl;
globalCellI = globalCells.toGlobal(localCellI);
}
@@ -130,6 +128,14 @@ const
<< exit(FatalError);
}
+
+ label procI = globalCells.whichProcID(globalCellI);
+ label procCellI = globalCells.toLocal(procI, globalCellI);
+
+ Info<< "Found point " << keepPoint << " in cell " << procCellI
+ << " on processor " << procI << endl;
+
+
if (globalCells.isLocal(globalCellI))
{
cellLabels[i] = localCellI;
diff --git a/src/meshTools/regionSplit/regionSplit.C b/src/meshTools/regionSplit/regionSplit.C
index f7815a16d81ab7bdd83cc4b0d54c86b6fa866cfd..2583c2ddfe8540b730b718c0ea17448301e9b98d 100644
--- a/src/meshTools/regionSplit/regionSplit.C
+++ b/src/meshTools/regionSplit/regionSplit.C
@@ -379,6 +379,7 @@ Foam::label Foam::regionSplit::calcLocalRegionSplit
Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
(
+ const bool doGlobalRegions,
const boolList& blockedFace,
const List<labelPair>& explicitConnections,
@@ -395,7 +396,7 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
cellRegion
);
- if (!Pstream::parRun())
+ if (!doGlobalRegions)
{
return autoPtr<globalIndex>(new globalIndex(nLocalRegions));
}
@@ -422,7 +423,7 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// (this will create gaps in the global region list so they will get
// merged later on)
- while (Pstream::parRun())
+ while (true)
{
if (debug)
{
@@ -690,13 +691,14 @@ Foam::autoPtr<Foam::globalIndex> Foam::regionSplit::calcRegionSplit
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
-Foam::regionSplit::regionSplit(const polyMesh& mesh)
+Foam::regionSplit::regionSplit(const polyMesh& mesh, const bool doGlobalRegions)
:
MeshObject<polyMesh, Foam::TopologicalMeshObject, regionSplit>(mesh),
labelList(mesh.nCells(), -1)
{
globalNumberingPtr_ = calcRegionSplit
(
+ doGlobalRegions, //do global regions
boolList(0, false), //blockedFaces
List<labelPair>(0), //explicitConnections,
*this
@@ -707,7 +709,8 @@ Foam::regionSplit::regionSplit(const polyMesh& mesh)
Foam::regionSplit::regionSplit
(
const polyMesh& mesh,
- const boolList& blockedFace
+ const boolList& blockedFace,
+ const bool doGlobalRegions
)
:
MeshObject<polyMesh, Foam::TopologicalMeshObject, regionSplit>(mesh),
@@ -715,6 +718,7 @@ Foam::regionSplit::regionSplit
{
globalNumberingPtr_ = calcRegionSplit
(
+ doGlobalRegions,
blockedFace, //blockedFaces
List<labelPair>(0), //explicitConnections,
*this
@@ -726,7 +730,8 @@ Foam::regionSplit::regionSplit
(
const polyMesh& mesh,
const boolList& blockedFace,
- const List<labelPair>& explicitConnections
+ const List<labelPair>& explicitConnections,
+ const bool doGlobalRegions
)
:
MeshObject<polyMesh, Foam::TopologicalMeshObject, regionSplit>(mesh),
@@ -734,6 +739,7 @@ Foam::regionSplit::regionSplit
{
globalNumberingPtr_ = calcRegionSplit
(
+ doGlobalRegions,
blockedFace, //blockedFaces
explicitConnections, //explicitConnections,
*this
diff --git a/src/meshTools/regionSplit/regionSplit.H b/src/meshTools/regionSplit/regionSplit.H
index b7be9dd32eb51567ad036050bb50f21d1f4adfee..83302529df24349890b4c4425d21fe4a0c5dfdae 100644
--- a/src/meshTools/regionSplit/regionSplit.H
+++ b/src/meshTools/regionSplit/regionSplit.H
@@ -87,6 +87,9 @@ Description
proc0 | proc1 | proc2
+ Can optionally keep all regions local to the processor.
+
+
SourceFiles
regionSplit.C
@@ -155,6 +158,7 @@ class regionSplit
//- Calculate global region split. Return globalIndex.
autoPtr<globalIndex> calcRegionSplit
(
+ const bool doGlobalRegions,
const boolList& blockedFace,
const List<labelPair>& explicitConnections,
labelList& cellRegion
@@ -170,11 +174,20 @@ public:
// Constructors
//- Construct from mesh
- regionSplit(const polyMesh&);
+ regionSplit
+ (
+ const polyMesh&,
+ const bool doGlobalRegions = Pstream::parRun()
+ );
//- Construct from mesh and whether face is blocked
// NOTE: blockedFace has to be consistent across coupled faces!
- regionSplit(const polyMesh&, const boolList& blockedFace);
+ regionSplit
+ (
+ const polyMesh&,
+ const boolList& blockedFace,
+ const bool doGlobalRegions = Pstream::parRun()
+ );
//- Construct from mesh and whether face is blocked. Additional explicit
// connections between normal boundary faces.
@@ -183,7 +196,8 @@ public:
(
const polyMesh&,
const boolList& blockedFace,
- const List<labelPair>&
+ const List<labelPair>&,
+ const bool doGlobalRegions = Pstream::parRun()
);
@@ -195,6 +209,12 @@ public:
return globalNumberingPtr_();
}
+ //- Return local number of regions
+ label nLocalRegions() const
+ {
+ return globalNumbering().localSize(Pstream::myProcNo());
+ }
+
//- Return total number of regions
label nRegions() const
{
diff --git a/src/meshTools/triSurface/booleanOps/surfaceIntersection/edgeIntersections.C b/src/meshTools/triSurface/booleanOps/surfaceIntersection/edgeIntersections.C
index a233af491d9e0f8b6e80a9d748891122d1843afc..441aa3be52b02497c4c97533b85f80af6bdbbdb4 100644
--- a/src/meshTools/triSurface/booleanOps/surfaceIntersection/edgeIntersections.C
+++ b/src/meshTools/triSurface/booleanOps/surfaceIntersection/edgeIntersections.C
@@ -117,19 +117,16 @@ void Foam::edgeIntersections::intersectEdges
<< " triangles ..." << endl;
}
- // Construct octree.
- const indexedOctree<treeDataTriSurface>& tree = querySurf2.tree();
-
-
- label nHits = 0;
-
+ pointField start(edgeLabels.size());
+ pointField end(edgeLabels.size());
+ vectorField edgeDirs(edgeLabels.size());
// Go through all edges, calculate intersections
forAll(edgeLabels, i)
{
label edgeI = edgeLabels[i];
- if (debug && (i % 1000 == 0))
+ if (debug)// && (i % 1000 == 0))
{
Pout<< "Intersecting edge " << edgeI << " with surface" << endl;
}
@@ -140,85 +137,81 @@ void Foam::edgeIntersections::intersectEdges
const point& pEnd = points1[meshPoints[e.end()]];
const vector eVec(pEnd - pStart);
- const scalar eMag = mag(eVec);
- const vector n(eVec/(eMag + VSMALL));
-
- // Smallish length for intersection calculation errors.
- const point tolVec = 1e-6*eVec;
+ const vector n(eVec/(mag(eVec) + VSMALL));
- // Start tracking somewhat before pStart and upto somewhat after p1.
+ // Start tracking somewhat before pStart and up to somewhat after p1.
// Note that tolerances here are smaller than those used to classify
// hit below.
// This will cause this hit to be marked as degenerate and resolved
// later on.
- point p0 = pStart - 0.5*surf1PointTol[e[0]]*n;
- const point p1 = pEnd + 0.5*surf1PointTol[e[1]]*n;
- const scalar maxS = mag(p1 - pStart);
+ start[i] = pStart - 0.5*surf1PointTol[e[0]]*n;
+ end[i] = pEnd + 0.5*surf1PointTol[e[1]]*n;
- // Get all intersections of the edge with the surface
+ edgeDirs[i] = n;
+ }
- DynamicList<pointIndexHit> currentIntersections(100);
- DynamicList<label> currentIntersectionTypes(100);
+ List<List<pointIndexHit> > edgeIntersections;
+ querySurf2.findLineAll
+ (
+ start,
+ end,
+ edgeIntersections
+ );
- while (true)
- {
- pointIndexHit pHit = tree.findLine(p0, p1);
+ label nHits = 0;
- if (pHit.hit())
- {
- nHits++;
+ // Classify the hits
+ forAll(edgeLabels, i)
+ {
+ const label edgeI = edgeLabels[i];
- currentIntersections.append(pHit);
+ labelList& intersectionTypes = classification_[edgeI];
+ intersectionTypes.setSize(edgeIntersections[i].size(), -1);
- // Classify point on surface1 edge.
- label edgeEnd = -1;
+ this->operator[](edgeI).transfer(edgeIntersections[i]);
- if (mag(pHit.hitPoint() - pStart) < surf1PointTol[e[0]])
- {
- edgeEnd = 0;
- }
- else if (mag(pHit.hitPoint() - pEnd) < surf1PointTol[e[1]])
- {
- edgeEnd = 1;
- }
- else if (mag(n & normals2[pHit.index()]) < alignedCos_)
- {
- Pout<< "Flat angle edge:" << edgeI
- << " face:" << pHit.index()
- << " cos:" << mag(n & normals2[pHit.index()])
- << endl;
- edgeEnd = 2;
- }
+ forAll(intersectionTypes, hitI)
+ {
+ const pointIndexHit& pHit = this->operator[](edgeI)[hitI];
+ label& hitType = intersectionTypes[hitI];
- currentIntersectionTypes.append(edgeEnd);
+ if (!pHit.hit())
+ {
+ continue;
+ }
- if (edgeEnd == 1)
- {
- // Close to end
- break;
- }
- else
- {
- // Continue tracking. Shift by a small amount.
- p0 = pHit.hitPoint() + tolVec;
+ const edge& e = surf1.edges()[edgeI];
- if (((p0-pStart) & n) >= maxS)
- {
- break;
- }
- }
+ // Classify point on surface1 edge.
+ if (mag(pHit.hitPoint() - start[i]) < surf1PointTol[e[0]])
+ {
+ // Intersection is close to edge start
+ hitType = 0;
}
- else
+ else if (mag(pHit.hitPoint() - end[i]) < surf1PointTol[e[1]])
{
- // No hit.
- break;
+ // Intersection is close to edge end
+ hitType = 1;
}
- }
+ else if (mag(edgeDirs[i] & normals2[pHit.index()]) < alignedCos_)
+ {
+ // Edge is almost coplanar with the face
+ Pout<< "Flat angle edge:" << edgeI
+ << " face:" << pHit.index()
+ << " cos:" << mag(edgeDirs[i] & normals2[pHit.index()])
+ << endl;
- // Done current edge. Transfer all data into *this
- operator[](edgeI).transfer(currentIntersections);
- classification_[edgeI].transfer(currentIntersectionTypes);
+ hitType = 2;
+ }
+
+ if (debug)
+ {
+ Info<< " hit " << pHit << " classify = " << hitType << endl;
+ }
+
+ nHits++;
+ }
}
if (debug)
@@ -226,7 +219,6 @@ void Foam::edgeIntersections::intersectEdges
Pout<< "Found " << nHits << " intersections of edges with surface ..."
<< endl;
}
-
}
diff --git a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
index b54ad7143d7adce09b0f47a6e392d347b3d054e4..0ede58107cbe86cc00e33457d6b73f8a5c8e1b64 100644
--- a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
+++ b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.C
@@ -202,7 +202,6 @@ void Foam::surfaceIntersection::storeIntersection
DynamicList<point>& allCutPoints
)
{
-
forAll(facesA, facesAI)
{
label faceA = facesA[facesAI];
@@ -969,10 +968,10 @@ Foam::surfaceIntersection::surfaceIntersection
if (!usedPoints.found(pointI))
{
- FatalErrorIn("surfaceIntersection::surfaceIntersection")
+ WarningIn("surfaceIntersection::surfaceIntersection")
<< "Problem: cut point:" << pointI
<< " coord:" << cutPoints_[pointI]
- << " not used by any edge" << abort(FatalError);
+ << " not used by any edge" << endl;
}
}
}
@@ -1147,6 +1146,12 @@ const Foam::edgeList& Foam::surfaceIntersection::cutEdges() const
}
+const Foam::labelPairLookup& Foam::surfaceIntersection::facePairToVertex() const
+{
+ return facePairToVertex_;
+}
+
+
const Foam::labelPairLookup& Foam::surfaceIntersection::facePairToEdge() const
{
return facePairToEdge_;
diff --git a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.H b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.H
index 68399c7e41d78d30501fb646b4054b131ba842c3..73638b43320f485072be761526fc50db108a5656 100644
--- a/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.H
+++ b/src/meshTools/triSurface/booleanOps/surfaceIntersection/surfaceIntersection.H
@@ -273,7 +273,7 @@ public:
const edgeList& cutEdges() const;
- //const labelPairLookup& facePairToVertex() const;
+ const labelPairLookup& facePairToVertex() const;
const labelPairLookup& facePairToEdge() const;
diff --git a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
index 80212847234373ee2b9e9254cd34f8e1e5b430f0..4e19610373bbfb9e5c9a0cd254767f244273f1ff 100644
--- a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
+++ b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.C
@@ -167,7 +167,7 @@ void Foam::decompositionMethod::calcCellCells
(
const polyMesh& mesh,
const labelList& agglom,
- const label nCoarse,
+ const label nLocalCoarse,
const bool parallel,
CompactListList<label>& cellCells
)
@@ -182,7 +182,7 @@ void Foam::decompositionMethod::calcCellCells
globalIndex globalAgglom
(
- nCoarse,
+ nLocalCoarse,
Pstream::msgType(),
Pstream::worldComm,
parallel
@@ -224,7 +224,7 @@ void Foam::decompositionMethod::calcCellCells
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// Number of faces per coarse cell
- labelList nFacesPerCell(nCoarse, 0);
+ labelList nFacesPerCell(nLocalCoarse, 0);
for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
{
@@ -374,7 +374,7 @@ void Foam::decompositionMethod::calcCellCells
// const boolList& blockedFace,
// const List<labelPair>& explicitConnections,
// const labelList& agglom,
-// const label nCoarse,
+// const label nLocalCoarse,
// const bool parallel,
// CompactListList<label>& cellCells
//)
@@ -389,7 +389,7 @@ void Foam::decompositionMethod::calcCellCells
//
// globalIndex globalAgglom
// (
-// nCoarse,
+// nLocalCoarse,
// Pstream::msgType(),
// Pstream::worldComm,
// parallel
@@ -431,7 +431,7 @@ void Foam::decompositionMethod::calcCellCells
// // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
//
// // Number of faces per coarse cell
-// labelList nFacesPerCell(nCoarse, 0);
+// labelList nFacesPerCell(nLocalCoarse, 0);
//
// // 1. Internal faces
// for (label faceI = 0; faceI < mesh.nInternalFaces(); faceI++)
@@ -764,6 +764,24 @@ Foam::labelList Foam::decompositionMethod::decompose
// Any weights specified?
label nWeights = returnReduce(cellWeights.size(), sumOp<label>());
+ if (nWeights > 0 && cellWeights.size() != mesh.nCells())
+ {
+ FatalErrorIn
+ (
+ "decompositionMethod::decompose\n"
+ "(\n"
+ " const polyMesh&,\n"
+ " const scalarField&,\n"
+ " const boolList&,\n"
+ " const PtrList<labelList>&,\n"
+ " const labelList&,\n"
+ " const List<labelPair>&\n"
+ ) << "Number of weights " << cellWeights.size()
+ << " differs from number of cells " << mesh.nCells()
+ << exit(FatalError);
+ }
+
+
// Any processor sets?
label nProcSets = 0;
forAll(specifiedProcessorFaces, setI)
@@ -828,9 +846,17 @@ Foam::labelList Foam::decompositionMethod::decompose
<< nProcSets << endl << endl;
}
- // Determine global regions, separated by blockedFaces
- regionSplit globalRegion(mesh, blockedFace, explicitConnections);
+ // Determine local regions, separated by blockedFaces
+ regionSplit localRegion(mesh, blockedFace, explicitConnections, false);
+
+ if (debug)
+ {
+ Info<< "Constrained decomposition:" << endl
+ << " split into " << localRegion.nLocalRegions()
+ << " regions."
+ << endl;
+ }
// Determine region cell centres
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
@@ -840,11 +866,11 @@ Foam::labelList Foam::decompositionMethod::decompose
// somewhere in the middle of the domain which might not be anywhere
// near any of the cells.
- pointField regionCentres(globalRegion.nRegions(), point::max);
+ pointField regionCentres(localRegion.nLocalRegions(), point::max);
- forAll(globalRegion, cellI)
+ forAll(localRegion, cellI)
{
- label regionI = globalRegion[cellI];
+ label regionI = localRegion[cellI];
if (regionCentres[regionI] == point::max)
{
@@ -855,22 +881,22 @@ Foam::labelList Foam::decompositionMethod::decompose
// Do decomposition on agglomeration
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
- scalarField regionWeights(globalRegion.nRegions(), 0);
+ scalarField regionWeights(localRegion.nLocalRegions(), 0);
if (nWeights > 0)
{
- forAll(globalRegion, cellI)
+ forAll(localRegion, cellI)
{
- label regionI = globalRegion[cellI];
+ label regionI = localRegion[cellI];
regionWeights[regionI] += cellWeights[cellI];
}
}
else
{
- forAll(globalRegion, cellI)
+ forAll(localRegion, cellI)
{
- label regionI = globalRegion[cellI];
+ label regionI = localRegion[cellI];
regionWeights[regionI] += 1.0;
}
@@ -879,7 +905,7 @@ Foam::labelList Foam::decompositionMethod::decompose
finalDecomp = decompose
(
mesh,
- globalRegion,
+ localRegion,
regionCentres,
regionWeights
);
diff --git a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.H b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.H
index 36bb6ab7a8a5a6ada21ab9be7265039928777a00..2d02936d2b1b9b3409cc4b606544e807bd07f338 100644
--- a/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.H
+++ b/src/parallel/decompose/decompositionMethods/decompositionMethod/decompositionMethod.H
@@ -219,7 +219,7 @@ public:
// Other
//- Helper: determine (local or global) cellCells from mesh
- // agglomeration.
+ // agglomeration. Agglomeration is local to the processor.
// local : connections are in local indices. Coupled across
// cyclics but not processor patches.
// global : connections are in global indices. Coupled across
@@ -228,34 +228,11 @@ public:
(
const polyMesh& mesh,
const labelList& agglom,
- const label nCoarse,
+ const label nLocalCoarse,
const bool global,
CompactListList<label>& cellCells
);
- //- Helper: determine (local or global) cellCells from mesh
- // agglomeration and additional specification:
- // - any additional connections between non-coupled internal
- // or boundary faces.
- // - any internal or coupled faces (or additional connections)
- // are blocked
- //
- // local : connections are in local indices. Coupled across
- // cyclics but not processor patches.
- // global : connections are in global indices. Coupled across
- // cyclics and processor patches.
- //static void calcCellCells
- //(
- // const polyMesh& mesh,
- // const boolList& blockedFace,
- // const List<labelPair>& explicitConnections,
- // const labelList& agglom,
- // const label nCoarse,
- // const bool global,
- // CompactListList<label>& cellCells
- //);
-
-
//- Helper: extract constraints:
// blockedface: existing faces where owner and neighbour on same
// proc
diff --git a/src/thermophysicalModels/laminarFlameSpeed/Gulders/Gulders.C b/src/thermophysicalModels/laminarFlameSpeed/Gulders/Gulders.C
index 9bdbaee7f65278826e1eb7ad2677ae51063aee96..051d1b2f1f873c901ba2e78ac80cf059df07cd73 100644
--- a/src/thermophysicalModels/laminarFlameSpeed/Gulders/Gulders.C
+++ b/src/thermophysicalModels/laminarFlameSpeed/Gulders/Gulders.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -217,7 +217,7 @@ Foam::laminarFlameSpeedModels::Gulders::operator()() const
dimensionedScalar
(
psiuReactionThermo_.lookup("stoichiometricAirFuelMassRatio")
- )*ft/((1 + SMALL) - ft)
+ )*ft/max(1 - ft, SMALL)
);
}
else
diff --git a/tutorials/mesh/foamyHexMesh/blob/Allrun b/tutorials/mesh/foamyHexMesh/blob/Allrun
index abcd8419bda49e95091ba65d24e96469a7645fef..3ef4e5b75e12b5550ea0f118190fb8c1473e727d 100755
--- a/tutorials/mesh/foamyHexMesh/blob/Allrun
+++ b/tutorials/mesh/foamyHexMesh/blob/Allrun
@@ -8,7 +8,13 @@ cd ${0%/*} || exit 1 # run from this directory
cp $FOAM_TUTORIALS/resources/geometry/blob.stl.gz constant/triSurface/
runApplication foamyHexMesh
+
runApplication collapseEdges -latestTime -collapseFaces
+mv log.collapseEdges log.collapseFaces
+
+runApplication collapseEdges -latestTime -collapseFaceSet indirectPatchFaces
+mv log.collapseEdges log.collapseFaceSet
+
runApplication checkMesh -latestTime -allGeometry -allTopology
diff --git a/tutorials/mesh/foamyHexMesh/blob/Allrun-parallel b/tutorials/mesh/foamyHexMesh/blob/Allrun-parallel
index 3113b19787164c7bc76f5c6d889b8eb45b971a5c..1674277920df26e29974bc492a4ac7c8e4dbdb76 100755
--- a/tutorials/mesh/foamyHexMesh/blob/Allrun-parallel
+++ b/tutorials/mesh/foamyHexMesh/blob/Allrun-parallel
@@ -14,7 +14,13 @@ runApplication blockMesh -region backgroundMeshDecomposition
runApplication decomposePar -region backgroundMeshDecomposition
runParallel foamyHexMesh $nProc
+
runParallel collapseEdges $nProc -latestTime -collapseFaces
+mv log.collapseEdges log.collapseFaces
+
+runParallel collapseEdges $nProc -latestTime -collapseFaceSet indirectPatchFaces
+mv log.collapseEdges log.collapseFaceSet
+
runParallel checkMesh $nProc -latestTime -allTopology -allGeometry
runApplication reconstructParMesh -latestTime
diff --git a/tutorials/mesh/foamyHexMesh/flange/Allrun b/tutorials/mesh/foamyHexMesh/flange/Allrun
index fe467dd3c693fd25392ff9d6023fd3f90bca3007..959eae6fd35131c8f94256d99102fc5a37861099 100755
--- a/tutorials/mesh/foamyHexMesh/flange/Allrun
+++ b/tutorials/mesh/foamyHexMesh/flange/Allrun
@@ -8,7 +8,13 @@ cd ${0%/*} || exit 1 # run from this directory
cp $FOAM_TUTORIALS/resources/geometry/flange.stl.gz constant/triSurface/
runApplication foamyHexMesh
+
runApplication collapseEdges -latestTime -collapseFaces
+mv log.collapseEdges log.collapseFaces
+
+runApplication collapseEdges -latestTime -collapseFaceSet indirectPatchFaces
+mv log.collapseEdges log.collapseFaceSet
+
runApplication checkMesh -latestTime -allGeometry -allTopology
diff --git a/tutorials/mesh/foamyHexMesh/flange/Allrun-parallel b/tutorials/mesh/foamyHexMesh/flange/Allrun-parallel
index 237296826d54b95fa01b1298bdc5425a333ca772..8f33af0d7e13c4926bb343ecdeee6cdfa3d49173 100755
--- a/tutorials/mesh/foamyHexMesh/flange/Allrun-parallel
+++ b/tutorials/mesh/foamyHexMesh/flange/Allrun-parallel
@@ -15,7 +15,13 @@ runApplication blockMesh -region backgroundMeshDecomposition
runApplication decomposePar -region backgroundMeshDecomposition
runParallel foamyHexMesh $nProc
+
runParallel collapseEdges $nProc -latestTime -collapseFaces
+mv log.collapseEdges log.collapseFaces
+
+runParallel collapseEdges $nProc -latestTime -collapseFaceSet indirectPatchFaces
+mv log.collapseEdges log.collapseFaceSet
+
runParallel checkMesh $nProc -latestTime -allTopology -allGeometry
runApplication reconstructParMesh -latestTime
diff --git a/tutorials/mesh/foamyHexMesh/mixerVessel/system/foamyHexMeshDict b/tutorials/mesh/foamyHexMesh/mixerVessel/system/foamyHexMeshDict
index 9ec059422315e72b2386b31245ede606f9a0b29e..dfbab846cb9ca0a233822d5c48fdce06f6a01df9 100644
--- a/tutorials/mesh/foamyHexMesh/mixerVessel/system/foamyHexMeshDict
+++ b/tutorials/mesh/foamyHexMesh/mixerVessel/system/foamyHexMeshDict
@@ -59,9 +59,9 @@ surfaceConformation
extendedFeatureEdgeMesh "vessel.extendedFeatureEdgeMesh";
regions
{
+ patch0 {}
patch1 {}
patch2 {}
- patch3 {}
}
}
@@ -118,7 +118,7 @@ motionControl
minimumCellSizeCoeff 0.1;
- maxRefinementIterations 0;
+ maxRefinementIterations 1;
maxSmoothingIterations 100;
@@ -142,7 +142,7 @@ motionControl
regions
{
- patch1
+ patch0
{
surfaceCellSizeFunction uniformValue;
uniformValueCoeffs
@@ -154,7 +154,7 @@ motionControl
uniformCoeffs{}
}
- patch2
+ patch1
{
priority 2;
surfaceCellSizeFunction uniformValue;
@@ -170,7 +170,7 @@ motionControl
}
}
- patch3
+ patch2
{
priority 2;
surfaceCellSizeFunction uniformValue;
diff --git a/tutorials/mesh/foamyHexMesh/simpleShapes/Allclean b/tutorials/mesh/foamyHexMesh/simpleShapes/Allclean
index 3e9b928240b8f003339d5ecd2d644f04a964560f..94024b2b8cef441e2ad9c7d66fe59a769e30855c 100755
--- a/tutorials/mesh/foamyHexMesh/simpleShapes/Allclean
+++ b/tutorials/mesh/foamyHexMesh/simpleShapes/Allclean
@@ -15,6 +15,8 @@ rm -r snapToSurface?.obj tetsToSnapTo.obj > /dev/null 2>&1
rm domain coneAndSphere > /dev/null 2>&1
+rm -rf 0/
+
cleanCase
# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/mesh/foamyHexMesh/simpleShapes/Allrun b/tutorials/mesh/foamyHexMesh/simpleShapes/Allrun
index bfc610ad80ccdf35c605f496710c98a6ce715870..497c84f91fb3677a96f2b50629b4b25b652be626 100755
--- a/tutorials/mesh/foamyHexMesh/simpleShapes/Allrun
+++ b/tutorials/mesh/foamyHexMesh/simpleShapes/Allrun
@@ -22,5 +22,7 @@ runApplication surfaceBooleanFeatures intersection \
runApplication foamyHexMesh
+runApplication collapseEdges -latestTime -collapseFaceSet indirectPatchFaces
+
# ----------------------------------------------------------------- end-of-file
diff --git a/wmake/rules/General/CGAL b/wmake/rules/General/CGAL
index 167a420ab6d938a0b176382123424ffc88b35576..ce87bd65a71a694782c6ea95dc203ac1678e014c 100644
--- a/wmake/rules/General/CGAL
+++ b/wmake/rules/General/CGAL
@@ -1,5 +1,4 @@
CGAL_INC = \
- -Wno-old-style-cast \
-I$(CGAL_ARCH_PATH)/include \
-I$(MPFR_ARCH_PATH)/include \
-I$(GMP_ARCH_PATH)/include \
diff --git a/wmake/rules/linux64Clang/c++ b/wmake/rules/linux64Clang/c++
index 1f78395afb28f385b735d8243eabb80b7de55bc0..1a0343b5cfc81819819cf26dd254bf26613998e1 100644
--- a/wmake/rules/linux64Clang/c++
+++ b/wmake/rules/linux64Clang/c++
@@ -3,6 +3,10 @@
# -Woverloaded-virtual may produce spurious warnings, disable for now
c++WARN = -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -Wno-overloaded-virtual -Wno-unused-comparison
+
+# Suppress CGAL warnings
+c++CGALWARN = -Wno-c++11-extensions -Wno-sometimes-uninitialized -Wno-mismatched-tags
+
CC = clang++ -m64
include $(RULES)/c++$(WM_COMPILE_OPTION)
diff --git a/wmake/rules/linuxClang/c++ b/wmake/rules/linuxClang/c++
index 7ae53634e5bc8d104c48ea90be3d4c5782e077a8..49849377a6214a1e2520318dcea97d2a01263abd 100644
--- a/wmake/rules/linuxClang/c++
+++ b/wmake/rules/linuxClang/c++
@@ -3,6 +3,10 @@
# -Woverloaded-virtual may produce spurious warnings, disable for now
c++WARN = -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -Wnon-virtual-dtor -Wno-overloaded-virtual -Wno-unused-comparison
+
+# Suppress CGAL warnings
+c++CGALWARN = -Wno-c++11-extensions -Wno-sometimes-uninitialized -Wno-mismatched-tags
+
CC = clang++ -m32
include $(RULES)/c++$(WM_COMPILE_OPTION)