1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
11 OpenFOAM is free software; you can redistribute it and/or modify it
12 under the terms of the GNU General Public License as published by the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
16 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
26 Reads surface and applies surface regioning to a mesh. Uses boundaryMesh
29 \*---------------------------------------------------------------------------*/
33 #include "boundaryMesh.H"
36 #include "polyTopoChange.H"
37 #include "polyModifyFace.H"
38 #include "globalMeshData.H"
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 // Adds empty patch if not yet there. Returns patchID.
45 label addPatch(polyMesh& mesh, const word& patchName)
47 label patchI = mesh.boundaryMesh().findPatchID(patchName);
51 const polyBoundaryMesh& patches = mesh.boundaryMesh();
53 List<polyPatch*> newPatches(patches.size() + 1);
57 // Copy all old patches
60 const polyPatch& pp = patches[i];
74 // Add zero-sized patch
85 mesh.removeBoundary();
86 mesh.addPatches(newPatches);
88 Pout<< "Created patch " << patchName << " at " << patchI << endl;
92 Pout<< "Reusing patch " << patchName << " at " << patchI << endl;
99 // Repatch single face. Return true if patch changed.
102 const polyMesh& mesh,
103 const boundaryMesh& bMesh,
104 const labelList& nearest,
105 const labelList& surfToMeshPatch,
107 polyTopoChange& meshMod
110 bool changed = false;
112 label bFaceI = faceI - mesh.nInternalFaces();
114 if (nearest[bFaceI] != -1)
116 // Use boundary mesh one.
117 label bMeshPatchID = bMesh.whichPatch(nearest[bFaceI]);
119 label patchID = surfToMeshPatch[bMeshPatchID];
121 if (patchID != mesh.boundaryMesh().whichPatch(faceI))
123 label own = mesh.faceOwner()[faceI];
125 label zoneID = mesh.faceZones().whichZone(faceI);
127 bool zoneFlip = false;
131 const faceZone& fZone = mesh.faceZones()[zoneID];
133 zoneFlip = fZone.flipMap()[fZone.whichFace(faceI)];
140 mesh.faces()[faceI],// modified face
141 faceI, // label of face being modified
145 patchID, // patch for face
146 false, // remove from zone
147 zoneID, // zone for face
148 zoneFlip // face flip in zone
165 int main(int argc, char *argv[])
167 argList::noParallel();
168 argList::validArgs.append("surface file");
169 argList::validOptions.insert("faceSet", "faceSet name");
170 argList::validOptions.insert("tol", "fraction of mesh size");
172 # include "setRootCase.H"
173 # include "createTime.H"
174 # include "createPolyMesh.H"
176 fileName surfName(args.additionalArgs()[0]);
178 Info<< "Reading surface from " << surfName << " ..." << endl;
180 bool readSet = args.options().found("faceSet");
185 setName = args.options()["faceSet"];
187 Info<< "Repatching only the faces in faceSet " << setName
188 << " according to nearest surface triangle ..." << endl;
192 Info<< "Patching all boundary faces according to nearest surface"
193 << " triangle ..." << endl;
196 scalar searchTol = 1E-3;
198 if (args.options().found("tol"))
200 searchTol = readScalar(IStringStream(args.options()["tol"])());
203 // Get search box. Anything not within this box will not be considered.
204 const boundBox& meshBb = mesh.globalData().bb();
206 const vector searchSpan(searchTol*(meshBb.max() - meshBb.min()));
208 Info<< "All boundary faces further away than " << searchTol
209 << " of mesh bounding box " << meshBb
210 << " will keep their patch label ..." << endl;
213 Info<< "Before patching:" << nl
214 << " patch\tsize" << endl;
216 forAll(mesh.boundaryMesh(), patchI)
218 Info<< " " << mesh.boundaryMesh()[patchI].name() << '\t'
219 << mesh.boundaryMesh()[patchI].size() << endl;
227 // Load in the surface.
228 bMesh.readTriSurface(surfName);
230 // Add all the boundaryMesh patches to the mesh.
231 const PtrList<boundaryPatch>& bPatches = bMesh.patches();
233 // Map from surface patch ( = boundaryMesh patch) to polyMesh patch
234 labelList patchMap(bPatches.size());
238 patchMap[i] = addPatch(mesh, bPatches[i].name());
241 // Obtain nearest face in bMesh for each boundary face in mesh that
242 // is within search span.
243 // Note: should only determine for faceSet if working with that.
244 labelList nearest(bMesh.getNearest(mesh, searchSpan));
247 // Dump unmatched faces to faceSet for debugging.
248 faceSet unmatchedFaces(mesh, "unmatchedFaces", nearest.size()/100);
250 forAll(nearest, bFaceI)
252 if (nearest[bFaceI] == -1)
254 unmatchedFaces.insert(mesh.nInternalFaces() + bFaceI);
258 Pout<< "Writing all " << unmatchedFaces.size()
259 << " unmatched faces to faceSet "
260 << unmatchedFaces.name()
263 unmatchedFaces.write();
267 polyTopoChange meshMod(mesh);
273 faceSet faceLabels(mesh, setName);
274 Info<< "Read " << faceLabels.size() << " faces to repatch ..." << endl;
276 forAllConstIter(faceSet, faceLabels, iter)
278 label faceI = iter.key();
280 if (repatchFace(mesh, bMesh, nearest, patchMap, faceI, meshMod))
288 forAll(nearest, bFaceI)
290 label faceI = mesh.nInternalFaces() + bFaceI;
292 if (repatchFace(mesh, bMesh, nearest, patchMap, faceI, meshMod))
299 Pout<< "Changed " << nChanged << " boundary faces." << nl << endl;
303 meshMod.changeMesh(mesh, false);
305 Info<< "After patching:" << nl
306 << " patch\tsize" << endl;
308 forAll(mesh.boundaryMesh(), patchI)
310 Info<< " " << mesh.boundaryMesh()[patchI].name() << '\t'
311 << mesh.boundaryMesh()[patchI].size() << endl;
318 // Write resulting mesh
319 Info << "Writing modified mesh to time " << runTime.value() << endl;
324 Info<< "End\n" << endl;
330 // ************************************************************************* //