diff --git a/applications/utilities/mesh/manipulation/createPatch/createPatch.C b/applications/utilities/mesh/manipulation/createPatch/createPatch.C
index e1e428471279431e6157f18f734aa4971bd36cef..f9cc4cce30752dbb52b5e94c2634d8f90051eb04 100644
--- a/applications/utilities/mesh/manipulation/createPatch/createPatch.C
+++ b/applications/utilities/mesh/manipulation/createPatch/createPatch.C
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
+ \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -26,8 +26,15 @@ Description
Utility to create patches out of selected boundary faces. Faces come either
from existing patches or from a faceSet.
+ More specifically it:
+ - creates new patches (from selected boundary faces). Synchronise faces
+ on coupled patches.
+ - synchronises points on coupled boundaries
+ - remove patches with 0 faces in them
+
\*---------------------------------------------------------------------------*/
+#include "syncTools.H"
#include "argList.H"
#include "polyMesh.H"
#include "Time.H"
@@ -48,6 +55,23 @@ namespace Foam
defineTemplateTypeNameAndDebug(IOPtrList<dictionary>, 0);
}
+// Combine operator to synchronise points. We choose point nearest to origin so
+// we can use e.g. great,great,great as null value.
+class nearestEqOp
+{
+
+public:
+
+ void operator()(vector& x, const vector& y) const
+ {
+ if (magSqr(y) < magSqr(x))
+ {
+ x = y;
+ }
+ }
+};
+
+
label getPatch(const polyBoundaryMesh& patches, const word& patchName)
{
@@ -55,7 +79,7 @@ label getPatch(const polyBoundaryMesh& patches, const word& patchName)
if (patchI == -1)
{
- FatalErrorIn("createPatch")
+ FatalErrorIn("createPatch(const polyBoundaryMesh&, const word&)")
<< "Cannot find source patch " << patchName
<< endl << "Valid patch names are " << patches.names()
<< exit(FatalError);
@@ -105,47 +129,174 @@ void changePatchID
// Filter out the empty patches.
void filterPatches(polyMesh& mesh)
{
- Info<< "Filtering empty patches." << endl;
-
const polyBoundaryMesh& patches = mesh.boundaryMesh();
- // Old patches.
+ // Patches to keep
DynamicList<polyPatch*> allPatches(patches.size());
- labelList compactPatchMap(patches.size());
+ label nOldPatches = returnReduce(patches.size(), sumOp<label>());
// Copy old patches.
forAll(patches, patchI)
{
const polyPatch& pp = patches[patchI];
- if (pp.size() > 0)
+ // Note: reduce possible since non-proc patches guaranteed in same order
+ if (!isA<processorPolyPatch>(pp))
{
- compactPatchMap[patchI] = allPatches.size();
+ if (returnReduce(pp.size(), sumOp<label>()) > 0)
+ {
+ allPatches.append
+ (
+ pp.clone
+ (
+ patches,
+ allPatches.size(),
+ pp.size(),
+ pp.start()
+ ).ptr()
+ );
+ }
+ else
+ {
+ Info<< "Removing empty patch " << pp.name() << " at position "
+ << patchI << endl;
+ }
+ }
+ }
+ // Copy non-empty processor patches
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
- allPatches.append
- (
- pp.clone
+ if (isA<processorPolyPatch>(pp))
+ {
+ if (pp.size() > 0)
+ {
+ allPatches.append
(
- patches,
- allPatches.size(),
- pp.size(),
- pp.start()
- ).ptr()
- );
+ pp.clone
+ (
+ patches,
+ allPatches.size(),
+ pp.size(),
+ pp.start()
+ ).ptr()
+ );
+ }
+ else
+ {
+ Info<< "Removing empty processor patch " << pp.name()
+ << " at position " << patchI << endl;
+ }
}
- else
+ }
+
+ label nAllPatches = returnReduce(allPatches.size(), sumOp<label>());
+ if (nAllPatches != nOldPatches)
+ {
+ Info<< "Removing patches." << endl;
+ allPatches.shrink();
+ mesh.removeBoundary();
+ mesh.addPatches(allPatches);
+ }
+ else
+ {
+ Info<< "No patches removed." << endl;
+ }
+}
+
+
+// Dump for all patches the current match
+void dumpCyclicMatch(const fileName& prefix, const polyMesh& mesh)
+{
+ const polyBoundaryMesh& patches = mesh.boundaryMesh();
+
+ forAll(patches, patchI)
+ {
+ if (isA<cyclicPolyPatch>(patches[patchI]))
{
- Info<< "Removing empty patch " << pp.name() << " at position "
- << patchI << endl;
+ const cyclicPolyPatch& cycPatch =
+ refCast<const cyclicPolyPatch>(patches[patchI]);
- compactPatchMap[patchI] = -1;
+ label halfSize = cycPatch.size()/2;
+
+ // Dump halves
+ {
+ OFstream str(prefix+cycPatch.name()+"_half0.obj");
+ meshTools::writeOBJ
+ (
+ str,
+ static_cast<faceList>
+ (
+ SubList<face>
+ (
+ cycPatch,
+ halfSize
+ )
+ ),
+ cycPatch.points()
+ );
+ }
+ {
+ OFstream str(prefix+cycPatch.name()+"_half1.obj");
+ meshTools::writeOBJ
+ (
+ str,
+ static_cast<faceList>
+ (
+ SubList<face>
+ (
+ cycPatch,
+ halfSize,
+ halfSize
+ )
+ ),
+ cycPatch.points()
+ );
+ }
+
+// cycPatch.writeOBJ
+// (
+// prefix+cycPatch.name()+"_half0.obj",
+// SubList<face>
+// (
+// cycPatch,
+// halfSize
+// ),
+// cycPatch.points()
+// );
+// cycPatch.writeOBJ
+// (
+// prefix+cycPatch.name()+"_half1.obj",
+// SubList<face>
+// (
+// cycPatch,
+// halfSize,
+// halfSize
+// ),
+// cycPatch.points()
+// );
+
+ // Lines between corresponding face centres
+ OFstream str(prefix+cycPatch.name()+"_match.obj");
+ label vertI = 0;
+
+ for (label faceI = 0; faceI < halfSize; faceI++)
+ {
+ const point& fc0 = mesh.faceCentres()[cycPatch.start()+faceI];
+ meshTools::writeOBJ(str, fc0);
+ vertI++;
+
+ label nbrFaceI = halfSize + faceI;
+ const point& fc1 = mesh.faceCentres()[cycPatch.start()+nbrFaceI];
+ meshTools::writeOBJ(str, fc1);
+ vertI++;
+
+ str<< "l " << vertI-1 << ' ' << vertI << nl;
+ }
}
}
- allPatches.shrink();
-
- mesh.removeBoundary();
- mesh.addPatches(allPatches);
}
@@ -153,99 +304,147 @@ void filterPatches(polyMesh& mesh)
int main(int argc, char *argv[])
{
- argList::noParallel();
argList::validOptions.insert("overwrite", "");
# include "setRootCase.H"
# include "createTime.H"
-# include "createPolyMesh.H"
- bool overwrite = args.options().found("overwrite");
+ const bool overwrite = args.options().found("overwrite");
- Info<< "Reading createPatchDict\n" << endl;
+ Info<< "Reading createPatchDict." << nl << endl;
- PtrList<dictionary> patchSources
+ IOdictionary dict
(
- IOdictionary
+ IOobject
(
- IOobject
- (
- "createPatchDict",
- runTime.system(),
- runTime,
- IOobject::MUST_READ,
- IOobject::NO_WRITE,
- false
- )
- ).lookup("patches")
+ "createPatchDict",
+ runTime.system(),
+ runTime,
+ IOobject::MUST_READ,
+ IOobject::NO_WRITE,
+ false
+ )
);
+
+ // Whether to synchronise points
+ const Switch pointSync(dict.lookup("pointSync"));
+
+
+ // Set the matching tolerance so we can read illegal meshes
+ scalar tol = readScalar(dict.lookup("matchTolerance"));
+ Info<< "Using relative tolerance " << tol
+ << " to match up faces and points" << nl << endl;
+ coupledPolyPatch::matchTol = tol;
+
+
+# include "createPolyMesh.H"
+
const polyBoundaryMesh& patches = mesh.boundaryMesh();
+ // If running parallel check same patches everywhere
+ patches.checkParallelSync(true);
- // 1. All all new patches
- // ~~~~~~~~~~~~~~~~~~~~~~
- // Old and new patches.
- DynamicList<polyPatch*> allPatches(patches.size() + patchSources.size());
+ dumpCyclicMatch("initial_", mesh);
- // Copy old patches.
- forAll(patches, patchI)
- {
- const polyPatch& pp = patches[patchI];
+ // Read patch construct info from dictionary
+ PtrList<dictionary> patchSources(dict.lookup("patches"));
- Info<< "Copying patch " << pp.name() << " at position "
- << patchI << endl;
- allPatches.append
- (
- pp.clone
- (
- patches,
- patchI,
- pp.size(),
- pp.start()
- ).ptr()
- );
- }
- forAll(patchSources, addedI)
+ // 1. Add all new patches
+ // ~~~~~~~~~~~~~~~~~~~~~~
+
+ if (patchSources.size() > 0)
{
- const dictionary& dict = patchSources[addedI];
+ // Old and new patches.
+ DynamicList<polyPatch*> allPatches(patches.size()+patchSources.size());
- word patchName(dict.lookup("name"));
+ label startFaceI = mesh.nInternalFaces();
- label destPatchI = patches.findPatchID(patchName);
+ // Copy old patches.
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
- word patchType(dict.lookup("type"));
+ if (!isA<processorPolyPatch>(pp))
+ {
+ allPatches.append
+ (
+ pp.clone
+ (
+ patches,
+ patchI,
+ pp.size(),
+ startFaceI
+ ).ptr()
+ );
+ startFaceI += pp.size();
+ }
+ }
- if (destPatchI == -1)
+ forAll(patchSources, addedI)
{
- destPatchI = allPatches.size();
+ const dictionary& dict = patchSources[addedI];
- Info<< "Adding new patch " << patchName << " of type " << patchType
- << " as patch " << destPatchI << endl;
+ word patchName(dict.lookup("name"));
- // Add an empty patch.
- allPatches.append
- (
- polyPatch::New
+ label destPatchI = patches.findPatchID(patchName);
+
+ word patchType(dict.lookup("type"));
+
+ if (destPatchI == -1)
+ {
+ destPatchI = allPatches.size();
+
+ Info<< "Adding new patch " << patchName
+ << " of type " << patchType
+ << " as patch " << destPatchI << endl;
+
+ // Add an empty patch.
+ allPatches.append
(
- patchType,
- patchName,
- 0,
- mesh.nFaces(),
- destPatchI,
- patches
- ).ptr()
- );
+ polyPatch::New
+ (
+ patchType,
+ patchName,
+ 0, // size
+ startFaceI, // start
+ destPatchI,
+ patches
+ ).ptr()
+ );
+ }
+ }
+
+ // Copy old patches.
+ forAll(patches, patchI)
+ {
+ const polyPatch& pp = patches[patchI];
+
+ if (isA<processorPolyPatch>(pp))
+ {
+ allPatches.append
+ (
+ pp.clone
+ (
+ patches,
+ patchI,
+ pp.size(),
+ startFaceI
+ ).ptr()
+ );
+ startFaceI += pp.size();
+ }
}
- }
- allPatches.shrink();
+ allPatches.shrink();
+ mesh.removeBoundary();
+ mesh.addPatches(allPatches);
- mesh.removeBoundary();
- mesh.addPatches(allPatches);
+ Info<< endl;
+ }
@@ -303,8 +502,8 @@ int main(int argc, char *argv[])
faceSet faces(mesh, setName);
- Info<< "Read " << faces.size() << " faces from faceSet "
- << faces.name() << endl;
+ Info<< "Read " << returnReduce(faces.size(), sumOp<label>())
+ << " faces from faceSet " << faces.name() << endl;
// Sort (since faceSet contains faces in arbitrary order)
labelList faceLabels(faces.toc());
@@ -341,15 +540,45 @@ int main(int argc, char *argv[])
<< "Valid source types are 'patches' 'set'" << exit(FatalError);
}
}
+ Info<< endl;
- // Change mesh, no inflation
- meshMod.changeMesh(mesh, false);
+ // Change mesh, use inflation to reforce calculation of transformation
+ // tensors.
+ Info<< "Doing topology modification to order faces." << nl << endl;
+ autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, true);
+ mesh.movePoints(map().preMotionPoints());
+ // Synchronise points.
+ if (!pointSync)
+ {
+ Info<< "Not synchronising points." << nl << endl;
+ }
+ else
+ {
+ Info<< "Synchronising points." << endl;
+
+ pointField newPoints(mesh.points());
+ syncTools::syncPointList
+ (
+ mesh,
+ newPoints,
+ nearestEqOp(), // cop
+ point(GREAT, GREAT, GREAT), // nullValue
+ true // applySeparation
+ );
+
+ scalarField diff(mag(newPoints-mesh.points()));
+ Info<< "Points changed by average:" << gAverage(diff)
+ << " max:" << gMax(diff) << nl << endl;
+
+ mesh.movePoints(newPoints);
+ }
// 3. Remove zeros-sized patches
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+ Info<< "Removing patches with no faces in them." << nl<< endl;
filterPatches(mesh);
@@ -363,7 +592,7 @@ int main(int argc, char *argv[])
}
// Write resulting mesh
- Info<< "Writing repatched mesh to " << runTime.timeName() << endl;
+ Info<< "Writing repatched mesh to " << runTime.timeName() << nl << endl;
mesh.write();
Info<< "End\n" << endl;
diff --git a/applications/utilities/mesh/manipulation/createPatch/createPatchDict b/applications/utilities/mesh/manipulation/createPatch/createPatchDict
index f11d33a0ea769e1af74ff70e9227beb1d7cccc4c..31d0078021782282424e8a4ff9b76556afcbb751 100644
--- a/applications/utilities/mesh/manipulation/createPatch/createPatchDict
+++ b/applications/utilities/mesh/manipulation/createPatch/createPatchDict
@@ -26,6 +26,12 @@ FoamFile
// face times this factor.
matchTolerance 1E-3;
+// Do a synchronisation of coupled points.
+pointSync true;
+
+
+// Patches to create.
+// If no patches does a coupled point and face synchronisation anyway.
patches
(
{