83 const label patchi =
bMesh.findPatchID(
name);
88 <<
"Known patches: " <<
bMesh.names() <<
endl;
93 if (
bMesh[patchi].empty())
95 Info<<
"Patch " <<
name <<
" has zero size" <<
nl;
104 int main(
int argc,
char *argv[])
108 "Merge the faces on specified patches (if geometrically possible)" 110 "faces become internal.\n" 111 "This utility can be called without arguments (uses stitchMeshDict)" 113 "two arguments (master/slave patch names)." 116 argList::noParallel();
117 argList::noFunctionObjects();
122 argList::addOption(
"dict",
"file",
"Alternative stitchMeshDict");
124 argList::addBoolOption
127 "Couple integral master/slave patches (2 argument mode: default)" 129 argList::addBoolOption
132 "Couple partially overlapping master/slave patches (2 argument mode)" 134 argList::addBoolOption
137 "Couple perfectly aligned master/slave patches (2 argument mode)" 139 argList::addBoolOption
142 "Write intermediate stages, not just the final result" 148 "Dictionary file with tolerances" 152 argList::noMandatoryArgs();
156 "The master patch name (non-dictionary mode)" 161 "The slave patch name (non-dictionary mode)" 167 const bool useCommandArgs = (
args.
size() > 1);
174 <<
"Cannot specify both dictionary and command-line arguments" 192 <<
"Only specify -integral with command-line arguments" 199 <<
"Only specify -partial with command-line arguments" 206 <<
"Only specify -perfect with command-line arguments" 214 const word oldInstance =
mesh.pointsInstance();
216 const bool intermediate =
args.
found(
"intermediate");
217 const bool overwrite =
args.
found(
"overwrite");
219 const word
dictName(
"stitchMeshDict");
227 const int integralCover =
args.
found(
"integral");
228 const int partialCover =
args.
found(
"partial");
229 const int perfectCover =
args.
found(
"perfect");
231 if ((integralCover + partialCover + perfectCover) > 1)
234 <<
"Can only specify one of -integral | -partial | -perfect." 236 <<
"Use perfect match option if the patches perfectly align" 237 <<
" (both vertex positions and face centres)" <<
endl 242 const word masterPatchName(
args[1]);
243 const word slavePatchName(
args[2]);
246 Info<<
" " << masterPatchName
247 <<
" / " << slavePatchName <<
nl;
268 dict.add(
"match", word(
"perfect"));
270 else if (partialCover)
275 slidingInterface::typeOfMatchNames[slidingInterface::PARTIAL]
283 slidingInterface::typeOfMatchNames[slidingInterface::INTEGRAL]
288 dict.add(
"master", masterPatchName);
289 dict.add(
"slave", slavePatchName);
291 validatedDict.add(
"stitchMesh",
dict);
301 IOdictionary stitchDict(
dictIO);
303 Info<<
" with " << stitchDict.size() <<
" entries" <<
nl;
308 for (
const entry& dEntry : stitchDict)
310 if (!dEntry.isDict())
312 Info<<
"Ignoring non-dictionary entry: " 313 << dEntry.keyword() <<
nl;
317 const keyType&
key = dEntry.keyword();
322 if (
dict.readIfPresent(
"match", matchName))
326 matchName !=
"perfect" 327 && !slidingInterface::typeOfMatchNames.
found(matchName)
330 Info<<
"Error: unknown match type - " << matchName
331 <<
" should be one of " 332 << slidingInterface::typeOfMatchNames.toc() <<
nl;
338 const word masterPatchName(
dict.get<word>(
"master"));
339 const word slavePatchName(
dict.get<word>(
"slave"));
342 Info<<
" " << masterPatchName
343 <<
" / " << slavePatchName <<
nl;
345 if (!requestedPatches.insert(masterPatchName))
347 Info<<
"Error: patch specified multiple times - " 348 << masterPatchName <<
nl;
352 if (!requestedPatches.insert(slavePatchName))
354 Info<<
"Error: patch specified multiple times - " 355 << slavePatchName <<
nl;
357 requestedPatches.erase(masterPatchName);
368 requestedPatches.erase(masterPatchName);
369 requestedPatches.erase(slavePatchName);
379 const label nActions = validatedDict.size();
381 Info<<
nl << nActions <<
" validated actions" <<
endl;
397 IOdictionary toleranceFile
401 args[
"toleranceDict"],
404 IOobject::MUST_READ_IF_MODIFIED,
408 slidingTolerances += toleranceFile;
416 Info<<
"Reading all current volfields" <<
endl;
417 PtrList<volScalarField> volScalarFields;
420 PtrList<volVectorField> volVectorFields;
423 PtrList<volSphericalTensorField> volSphericalTensorFields;
426 PtrList<volSymmTensorField> volSymmTensorFields;
429 PtrList<volTensorField> volTensorFields;
444 IOstream::defaultPrecision(
max(10u, IOstream::defaultPrecision()));
446 polyTopoChanger stitcher(
mesh, IOobject::NO_READ);
450 for (
const entry& dEntry : validatedDict)
455 bool perfect =
false;
456 slidingInterface::typeOfMatch matchType = slidingInterface::PARTIAL;
459 if (
dict.readIfPresent(
"match", matchName))
461 if (matchName ==
"perfect")
467 matchType = slidingInterface::typeOfMatchNames[matchName];
472 const word masterPatchName(
dict.get<word>(
"master"));
473 const word slavePatchName(
dict.get<word>(
"slave"));
476 const word mergePatchName(masterPatchName + slavePatchName);
477 const word cutZoneName(mergePatchName +
"CutFaceZone");
479 Info<<
nl <<
"========================================" <<
nl;
484 Info<<
"Coupling PERFECTLY aligned patches " 485 << masterPatchName <<
" / " << slavePatchName <<
nl <<
nl 486 <<
"Resulting (internal) faces in faceZone " 487 << cutZoneName <<
nl <<
nl 488 <<
"The patch vertices and face centres must align within a" 489 <<
" tolerance relative to the minimum edge length on the patch" 492 else if (matchType == slidingInterface::INTEGRAL)
494 Info<<
"Coupling INTEGRALLY matching of patches " 495 << masterPatchName <<
" / " << slavePatchName <<
nl <<
nl 496 <<
"Resulting (internal) faces in faceZone " 497 << cutZoneName <<
nl <<
nl 498 <<
"The overall area covered by both patches should be" 499 <<
" identical!" <<
endl 500 <<
"If this is not the case use partial" 505 Info<<
"Coupling PARTIALLY overlapping patches " 506 << masterPatchName <<
" / " << slavePatchName <<
nl 507 <<
"Resulting internal faces in faceZone " 509 <<
"Uncovered faces remain in their patch" 515 const polyPatch& masterPatch =
mesh.boundaryMesh()[masterPatchName];
516 const polyPatch& slavePatch =
mesh.boundaryMesh()[slavePatchName];
518 mesh.pointZones().clearAddressing();
519 mesh.faceZones().clearAddressing();
520 mesh.cellZones().clearAddressing();
539 ).resetAddressing(std::move(faceIds),
false);
561 mergePatchName +
"CutPointZone",
567 mergePatchName +
"Side0Zone",
569 ).resetAddressing(std::move(faceIds),
false);
572 faceIds.resize_nocopy(slavePatch.size());
577 mergePatchName +
"Side1Zone",
579 ).resetAddressing(std::move(faceIds),
false);
598 mergePatchName +
"Side0Zone",
599 mergePatchName +
"Side1Zone",
600 mergePatchName +
"CutPointZone",
609 static_cast<slidingInterface&
>(stitcher[0]).setTolerances
619 if (!overwrite && (intermediate || actioni == nActions))
625 autoPtr<mapPolyMesh> morphMap = stitcher.changeMesh(
true);
627 mesh.movePoints(morphMap->preMotionPoints());
632 mesh.setInstance(oldInstance);
633 stitcher.instance() = oldInstance;
636 if (intermediate || actioni == nActions)
638 Info<<
nl <<
"Writing polyMesh to time " 644 !
runTime.objectRegistry::writeObject
656 <<
"Failed writing polyMesh." 660 mesh.faceZones().write();
661 mesh.pointZones().write();
662 mesh.cellZones().write();
wordList ReadFields(const typename GeoMesh::Mesh &mesh, const IOobjectList &objects, PtrList< GeometricField< Type, PatchField, GeoMesh >> &fields, const bool syncPar=true, const bool readOldTime=false)
Read Geometric fields of templated type.
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
constexpr char nl
The newline '\n' character (0x0a)
const word dictName("faMeshDefinition")
Ostream & endl(Ostream &os)
Add newline and flush stream.
void resize_nocopy(const label len)
Adjust allocated size of list without necessarily.
Field reading functions for post-processing utilities.
void exit(const int errNo=1)
Exit : can be called for any error to exit program.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
void clear()
Clear the list, i.e. set size to zero.
label size() const noexcept
The number of arguments.
HashSet< word, Hash< word > > wordHashSet
A HashSet of words, uses string hasher.
Ostream & flush(Ostream &os)
Flush stream.
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
messageStream Info
Information stream (stdout output on master, null elsewhere)
List< label > labelList
A List of labels.
Foam::argList args(argc, argv)
PrimitivePatch< List< face >, const pointField > bMesh
Holder of faceList and points. (v.s. e.g. primitivePatch which references points) ...
bool checkPatch(const bool allGeometry, const std::string &name, const polyMesh &mesh, const PatchType &pp, const labelUList &meshEdges, labelHashSet *pointSetPtr=nullptr, labelHashSet *badEdgesPtr=nullptr)
bool check(bool checkArgs=argList::argsMandatory(), bool checkOpts=true) const
Check the parsed command-line for mandatory arguments and that all the options are correct...
bool found(const word &optName) const
Return true if the named option is found.