113 bool changed =
false;
117 if (nearest[bFacei] != -1)
122 label patchID = surfToMeshPatch[rMeshPatchID];
130 bool zoneFlip =
false;
167 int main(
int argc,
char *argv[])
171 "reads surface and applies surface regioning to a mesh"
180 "only repatch the faces in specified faceSet"
186 "search tolerance as fraction of mesh size (default 1e-3)"
195 Info<<
"Reading surface from " << surfName <<
" ..." <<
endl;
202 Info<<
"Repatching only the faces in faceSet " << setName
203 <<
" according to nearest surface triangle ..." <<
endl;
207 Info<<
"Patching all boundary faces according to nearest surface"
208 <<
" triangle ..." <<
endl;
217 Info<<
"All boundary faces further away than " << searchTol
218 <<
" of mesh bounding box " <<
meshBb
219 <<
" will keep their patch label ..." <<
endl;
222 Info<<
"Before patching:" <<
nl
223 <<
" patch\tsize" <<
endl;
246 patchMap[i] = addPatch(mesh, bPatches[i].
name());
256 faceSet unmatchedFaces(mesh,
"unmatchedFaces", nearest.
size()/100);
260 if (nearest[bFacei] == -1)
266 Pout<<
"Writing all " << unmatchedFaces.size()
267 <<
" unmatched faces to faceSet "
268 << unmatchedFaces.
name()
271 unmatchedFaces.
write();
281 faceSet faceLabels(mesh, setName);
282 Info<<
"Read " << faceLabels.size() <<
" faces to repatch ..." <<
endl;
286 label facei = iter.key();
288 if (repatchFace(mesh, rMesh, nearest, patchMap, facei, meshMod))
300 if (repatchFace(mesh, rMesh, nearest, patchMap, facei, meshMod))
307 Pout<<
"Changed " << nChanged <<
" boundary faces." <<
nl <<
endl;
313 Info<<
"After patching:" <<
nl
314 <<
" patch\tsize" <<
endl;
327 Info<<
"Writing modified mesh to time "
328 << runTime.userTimeName() <<
endl;
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
const word & name() const
Return name.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
void size(const label)
Override size to be inconsistent with allocated storage.
label whichZone(const label objectIndex) const
Given a global object index, return the zone it is in.
virtual const fileName & name() const
Return the name of the stream.
virtual Ostream & write(const char)=0
Write character.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
label size() const
Return the number of elements in the UPtrList.
static void addOption(const word &opt, const string ¶m="", const string &usage="")
Add to an option to validOptions with usage information.
static void addNote(const string &)
Add extra notes for the usage information.
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
static void noParallel()
Remove the parallel options.
static SLList< string > validArgs
A list of valid (mandatory) arguments.
T optionLookupOrDefault(const word &opt, const T &deflt) const
Read a value from the named option if present.
A bounding box defined in terms of the points at its extremities.
A subset of mesh faces organised as a primitive patch.
const boolList & flipMap() const
Return face flip map.
label whichFace(const label globalCellID) const
Helper function to re-direct to zone::localID(...)
A class for handling file names.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
label findPatchID(const word &patchName) const
Find patch index given a name.
Mesh consisting of general polyhedral cells.
const meshFaceZones & faceZones() const
Return face zones.
virtual const faceList & faces() const
Return raw faces.
virtual const labelList & faceOwner() const
Return face owner.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
void addPatches(const List< polyPatch * > &, const bool validBoundary=true)
Add boundary patches.
void removeBoundary()
Remove boundary patches.
const boundBox & bounds() const
Return mesh bounding box.
Class describing modification of a face.
A patch is a list of labels that address the faces in the global face list.
label start() const
Return start label of this patch in the polyMesh face list.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
Direct mesh changes based on v1.3 polyTopoChange syntax.
autoPtr< polyTopoChangeMap > changeMesh(polyMesh &mesh, const bool inflate, const bool syncParallel=true, const bool orderCells=false, const bool orderPoints=false)
Inplace changes mesh without change of patches.
label setAction(const topoAction &action)
For compatibility with polyTopoChange: set topological action.
label nInternalFaces() const
virtual bool write(const bool write=true) const
Write using setting from DB.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
void readTriSurface(const fileName &)
Read from triSurface.
labelList getNearest(const primitiveMesh &pMesh, const vector &searchSpan) const
Get rMesh index of nearest face for every boundary face in.
label whichPatch(const label facei) const
Get index of patch face is in.
const PtrList< repatchPatch > & patches() const
Access the patches.
A class for handling words, derived from string.
int main(int argc, char *argv[])
const fvPatchList & patches
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
prefixOSstream Pout(cout, "Pout")
List< treeBoundBox > meshBb(1, treeBoundBox(boundBox(coarseMesh.points(), false)).extend(1e-3))
Foam::argList args(argc, argv)