surfaceToPatch.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2020-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Application
28 surfaceToPatch
29
30Group
31 grpSurfaceUtilities
32
33Description
34 Reads surface and applies surface regioning to a mesh. Uses boundaryMesh
35 to do the hard work.
36
37\*---------------------------------------------------------------------------*/
38
39#include "argList.H"
40#include "Time.H"
41#include "boundaryMesh.H"
42#include "polyMesh.H"
43#include "faceSet.H"
44#include "polyTopoChange.H"
45#include "polyModifyFace.H"
46#include "globalMeshData.H"
47
48using namespace Foam;
49
50// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51
52// Adds empty patch if not yet there. Returns patchID.
53label addPatch(polyMesh& mesh, const word& patchName)
54{
55 label patchi = mesh.boundaryMesh().findPatchID(patchName);
56
57 if (patchi == -1)
58 {
59 const polyBoundaryMesh& patches = mesh.boundaryMesh();
60
61 List<polyPatch*> newPatches(patches.size() + 1);
62
63 patchi = 0;
64
65 // Copy all old patches
66 forAll(patches, i)
67 {
68 const polyPatch& pp = patches[i];
69
70 newPatches[patchi] =
71 pp.clone
72 (
73 patches,
74 patchi,
75 pp.size(),
76 pp.start()
77 ).ptr();
78
79 patchi++;
80 }
81
82 // Add zero-sized patch
83 newPatches[patchi] =
84 new polyPatch
85 (
86 patchName,
87 0,
88 mesh.nFaces(),
89 patchi,
90 patches,
91 polyPatch::typeName
92 );
93
94 mesh.removeBoundary();
95 mesh.addPatches(newPatches);
96
97 Pout<< "Created patch " << patchName << " at " << patchi << endl;
98 }
99 else
100 {
101 Pout<< "Reusing patch " << patchName << " at " << patchi << endl;
102 }
103
104 return patchi;
105}
106
107
108// Repatch single face. Return true if patch changed.
109bool repatchFace
110(
111 const polyMesh& mesh,
112 const boundaryMesh& bMesh,
113 const labelList& nearest,
114 const labelList& surfToMeshPatch,
115 const label facei,
116 polyTopoChange& meshMod
117)
118{
119 bool changed = false;
120
121 label bFacei = facei - mesh.nInternalFaces();
122
123 if (nearest[bFacei] != -1)
124 {
125 // Use boundary mesh one.
126 label bMeshPatchID = bMesh.whichPatch(nearest[bFacei]);
127
128 label patchID = surfToMeshPatch[bMeshPatchID];
129
130 if (patchID != mesh.boundaryMesh().whichPatch(facei))
131 {
132 label own = mesh.faceOwner()[facei];
133
134 label zoneID = mesh.faceZones().whichZone(facei);
135
136 bool zoneFlip = false;
137
138 if (zoneID >= 0)
139 {
140 const faceZone& fZone = mesh.faceZones()[zoneID];
141
142 zoneFlip = fZone.flipMap()[fZone.whichFace(facei)];
143 }
144
145 meshMod.setAction
146 (
148 (
149 mesh.faces()[facei],// modified face
150 facei, // label of face being modified
151 own, // owner
152 -1, // neighbour
153 false, // face flip
154 patchID, // patch for face
155 false, // remove from zone
156 zoneID, // zone for face
157 zoneFlip // face flip in zone
158 )
159 );
160
161 changed = true;
162 }
163 }
164 else
165 {
166 changed = false;
167 }
168 return changed;
169}
170
171
172
173int main(int argc, char *argv[])
174{
175 argList::addNote
176 (
177 "Reads surface and applies surface regioning to a mesh"
178 );
179
180 argList::noParallel();
181 argList::addArgument("surfaceFile");
182 argList::addOption
183 (
184 "faceSet",
185 "name",
186 "Only repatch the faces in specified faceSet"
187 );
188 argList::addOption
189 (
190 "tol",
191 "scalar",
192 "Search tolerance as fraction of mesh size (default 1e-3)"
193 );
194
195 #include "setRootCase.H"
196 #include "createTime.H"
197 #include "createPolyMesh.H"
198
199 const auto surfName = args.get<fileName>(1);
200
201 Info<< "Reading surface from " << surfName << " ..." << endl;
202
203 word setName;
204 const bool readSet = args.readIfPresent("faceSet", setName);
205
206 if (readSet)
207 {
208 Info<< "Repatching only the faces in faceSet " << setName
209 << " according to nearest surface triangle ..." << endl;
210 }
211 else
212 {
213 Info<< "Patching all boundary faces according to nearest surface"
214 << " triangle ..." << endl;
215 }
216
217 const scalar searchTol = args.getOrDefault<scalar>("tol", 1e-3);
218
219 // Get search box. Anything not within this box will not be considered.
220 const boundBox& meshBb = mesh.bounds();
221 const vector searchSpan = searchTol * meshBb.span();
222
223 Info<< "All boundary faces further away than " << searchTol
224 << " of mesh bounding box " << meshBb
225 << " will keep their patch label ..." << endl;
226
227
228 Info<< "Before patching:" << nl
229 << " patch\tsize" << endl;
230
231 forAll(mesh.boundaryMesh(), patchi)
232 {
233 Info<< " " << mesh.boundaryMesh()[patchi].name() << '\t'
234 << mesh.boundaryMesh()[patchi].size() << nl;
235 }
236 Info<< endl;
237
238
240
241 // Load in the surface.
242 bMesh.readTriSurface(surfName);
243
244 // Add all the boundaryMesh patches to the mesh.
245 const PtrList<boundaryPatch>& bPatches = bMesh.patches();
246
247 // Map from surface patch ( = boundaryMesh patch) to polyMesh patch
248 labelList patchMap(bPatches.size());
249
250 forAll(bPatches, i)
251 {
252 patchMap[i] = addPatch(mesh, bPatches[i].name());
253 }
254
255 // Obtain nearest face in bMesh for each boundary face in mesh that
256 // is within search span.
257 // Note: should only determine for faceSet if working with that.
258 labelList nearest(bMesh.getNearest(mesh, searchSpan));
259
260 {
261 // Dump unmatched faces to faceSet for debugging.
262 faceSet unmatchedFaces(mesh, "unmatchedFaces", nearest.size()/100);
263
264 forAll(nearest, bFacei)
265 {
266 if (nearest[bFacei] == -1)
267 {
268 unmatchedFaces.insert(mesh.nInternalFaces() + bFacei);
269 }
270 }
271
272 Pout<< "Writing all " << unmatchedFaces.size()
273 << " unmatched faces to faceSet "
274 << unmatchedFaces.name()
275 << endl;
276
277 unmatchedFaces.write();
278 }
279
280
281 polyTopoChange meshMod(mesh);
282
283 label nChanged = 0;
284
285 if (readSet)
286 {
287 faceSet faceLabels(mesh, setName);
288 Info<< "Read " << faceLabels.size() << " faces to repatch ..." << endl;
289
290 for (const label facei : faceLabels)
291 {
292 if (repatchFace(mesh, bMesh, nearest, patchMap, facei, meshMod))
293 {
294 nChanged++;
295 }
296 }
297 }
298 else
299 {
300 forAll(nearest, bFacei)
301 {
302 const label facei = mesh.nInternalFaces() + bFacei;
303
304 if (repatchFace(mesh, bMesh, nearest, patchMap, facei, meshMod))
305 {
306 ++nChanged;
307 }
308 }
309 }
310
311 Pout<< "Changed " << nChanged << " boundary faces." << nl << endl;
312
313 if (nChanged > 0)
314 {
315 meshMod.changeMesh(mesh, false);
316
317 Info<< "After patching:" << nl
318 << " patch\tsize" << endl;
319
320 forAll(mesh.boundaryMesh(), patchi)
321 {
322 Info<< " " << mesh.boundaryMesh()[patchi].name() << '\t'
323 << mesh.boundaryMesh()[patchi].size() << endl;
324 }
325 Info<< endl;
326
327
328 ++runTime;
329
330 // Write resulting mesh
331 Info<< "Writing modified mesh to time " << runTime.value() << endl;
332 mesh.write();
333 }
334
335
336 Info<< "End\n" << endl;
337
338 return 0;
339}
340
341
342// ************************************************************************* //
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
virtual const fileName & name() const
Get the name of the stream.
Definition: OSstream.H:107
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
A list of faces which address into the list of points.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
label size() const noexcept
The number of elements in the list.
Definition: UPtrListI.H:106
T get(const label index) const
Get a value from the argument at index.
Definition: argListI.H:278
bool readIfPresent(const word &optName, T &val) const
Read a value from the named option if present.
Definition: argListI.H:323
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Definition: argListI.H:307
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:64
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
Definition: boundaryMesh.H:63
A list of face labels.
Definition: faceSet.H:54
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:67
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
Definition: faceZone.C:372
const boolList & flipMap() const noexcept
Return face flip map.
Definition: faceZone.H:272
A class for handling file names.
Definition: fileName.H:76
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Class describing modification of a face.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
virtual autoPtr< polyPatch > clone(const labelList &faceCells) const
Construct and return a clone, setting faceCells.
Definition: polyPatch.H:238
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:364
Direct mesh changes based on v1.3 polyTopoChange syntax.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
autoPtr< mapPolyMesh > changeMesh(polyMesh &mesh, const labelUList &patchMap, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
Definition: word.H:68
const polyBoundaryMesh & patches
dynamicFvMesh & mesh
engineTime & runTime
Required Variables.
const labelIOList & zoneID
Namespace for OpenFOAM.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
List< treeBoundBox > meshBb(1, treeBoundBox(boundBox(coarseMesh.points(), false)).extend(rndGen, 1e-3))
Foam::argList args(argc, argv)
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333