83 void createDummyFvMeshFiles(
const polyMesh& mesh,
const word& regionName)
98 Info<<
"Testing:" << io.objectPath() <<
endl;
102 Info<<
"Writing dummy " << regionName/io.name() <<
endl;
105 dummyDict.
add(
"divSchemes", divDict);
107 dummyDict.
add(
"gradSchemes", gradDict);
109 dummyDict.
add(
"laplacianSchemes", laplDict);
118 mesh.time().system(),
128 Info<<
"Writing dummy " << regionName/io.name() <<
endl;
143 <<
"Cannot find patch " << name
144 <<
" in the source mesh.\n" 145 <<
"Valid patch names are " << patches.
names()
158 const polyPatch& pp = patches[findPatchID(patches, names[i])];
166 const polyPatch& pp = patches[findPatchID(patches, names[i])];
170 faceLabels[n++] = pp.
start()+j;
187 label oldFacei = faceLabels[i];
189 if (reverseMap[oldFacei] >= 0)
191 newFaceLabels[newI++] = reverseMap[oldFacei];
207 label oldCelli = cellLabels[i];
209 if (reverseMap[oldCelli] >= 0)
211 newCellLabels.
insert(reverseMap[oldCelli]);
218 template<
class PatchType>
219 void changeFrontBackPatches
222 const word& frontPatchName,
223 const word& backPatchName
228 label frontPatchi = findPatchID(patches, frontPatchName);
229 label backPatchi = findPatchID(patches, backPatchName);
237 if (patchi == frontPatchi || patchi == backPatchi)
264 int main(
int argc,
char *argv[])
278 Info<<
"Create mesh " << regionName <<
" for time = " 279 << runTimeExtruded.timeName() <<
nl <<
endl;
284 Info<<
"Create mesh for time = " 285 << runTimeExtruded.timeName() <<
nl <<
endl;
297 const ExtrudeMode mode = ExtrudeModeNames.
read 305 Info<<
"Extruding from " << ExtrudeModeNames[
mode]
306 <<
" using model " << model().type() <<
endl;
309 Info<<
"Flipping normals before extruding" <<
endl;
313 Info<<
"Collapsing edges < " << mergeTol <<
" of bounding box" <<
endl;
317 Info<<
"Not collapsing any edges after extrusion" <<
endl;
335 if (mode == PATCH || mode == MESH)
337 if (flipNormals && mode == MESH)
340 <<
"Flipping normals not supported for extrusions from mesh." 355 dict.
lookup(
"sourcePatches") >> sourcePatches;
357 if (sourcePatches.
size() == 1)
359 frontPatchName = sourcePatches[0];
362 Info<<
"Extruding patches " << sourcePatches
363 <<
" on mesh " << sourceCasePath <<
nl 382 const labelList meshFaces(patchFaces(patches, sourcePatches));
384 if (mode == PATCH && flipNormals)
391 label meshFacei = meshFaces[i];
402 bool zoneFlip =
false;
406 zoneFlip = mesh.
faceZones()[zoneI].flipMap()[index];
411 mesh.
faces()[meshFacei].reverseFace(),
429 if (map().hasMotionPoints())
457 extrudePatch.meshEdges
512 Info<<
"Adding overall " << nAdded <<
" processor patches." <<
endl;
536 Pout<<
"Adding patch " << patchi
538 <<
" and " << nbrProci
556 mesh.removeFvBoundary();
557 mesh.addFvPatches(newPatches,
true);
566 dict.
lookup(
"exposedPatchName") >> backPatchName;
570 findPatchID(patches, backPatchName)
575 pointField layer0Points(extrudePatch.nPoints());
576 pointField displacement(extrudePatch.nPoints());
577 forAll(displacement, pointi)
579 const vector& patchNormal = extrudePatchPointNormals[pointi];
582 layer0Points[pointi] = model()
584 extrudePatch.localPoints()[pointi],
589 point extrudePt = model()
591 extrudePatch.localPoints()[pointi],
595 displacement[pointi] = extrudePt - layer0Points[pointi];
601 if (
gMax(
mag(layer0Points-extrudePatch.localPoints())) > small)
603 Info<<
"Moving mesh to layer0 points since differ from original" 604 <<
" points - this can happen for wedge extrusions." <<
nl 608 forAll(extrudePatch.meshPoints(), i)
610 newPoints[extrudePatch.meshPoints()[i]] = layer0Points[i];
617 labelList nFaceLayers(extrudePatch.size(), model().nLayers());
620 labelList nPointLayers(extrudePatch.nPoints(), model().nLayers());
623 vectorField firstLayerDisp(displacement*model().sumThickness(1));
639 layerExtrude.setRefinement
655 forAll(layerExtrude.addedPoints(), pointi)
657 const labelList& pPoints = layerExtrude.addedPoints()[pointi];
660 label meshPointi = pPoints[pPointi];
666 extrudePatch.localPoints()[pointi],
667 extrudePatchPointNormals[pointi],
675 )[meshPointi] = modelPt;
684 forAll(backPatchFaces, patchFacei)
686 backPatchFaces[patchFacei] = layerFaces[patchFacei].
first();
687 frontPatchFaces[patchFacei] = layerFaces[patchFacei].
last();
692 createDummyFvMeshFiles(mesh, regionDir);
701 runTimeExtruded.constant(),
710 layerExtrude.updateMesh
720 map().reverseFaceMap(),
725 map().reverseFaceMap(),
734 layerExtrude.addedCells
737 layerExtrude.layerFaces()
742 const labelList& aCells = addedCells[facei];
745 addedCellsSet.
insert(aCells[i]);
756 Info<<
"Extruding surfaceMesh read from file " << surfName <<
nl 767 faces[i] = fMesh[i].reverseFace();
771 Info<<
"Extruding surface with :" <<
nl 772 <<
" points : " << fMesh.points().size() <<
nl 773 <<
" faces : " << fMesh.size() <<
nl 774 <<
" normals[0] : " << fMesh.faceNormals()[0]
778 meshFromSurface.
reset 785 runTimeExtruded.constant(),
795 frontPatchName =
"originalPatch";
796 frontPatchFaces = patchFaces
798 meshFromSurface().boundaryMesh(),
801 backPatchName =
"otherSide";
802 backPatchFaces = patchFaces
804 meshFromSurface().boundaryMesh(),
820 const scalar mergeDim = mergeTol * bb.
minDim();
822 Info<<
"Mesh bounding box : " << bb <<
nl 823 <<
" with span : " << span <<
nl 824 <<
"Merge distance : " << mergeDim <<
nl 833 Pout<<
"Collapsing edges < " << mergeDim <<
" ..." <<
nl <<
endl;
846 const edge& e = edges[edgeI];
848 scalar d = e.
mag(points);
852 Pout<<
"Merging edge " << e <<
" since length " << d
853 <<
" << " << mergeDim <<
nl;
856 collapsePointToLocation.set(e[1], points[e[0]]);
864 collapser.consistentCollapse
868 collapsePointToLocation,
877 bool anyChange = collapser.setRefinement(allPointInfo, meshMod);
889 updateFaceLabels(map(), frontPatchFaces);
890 updateFaceLabels(map(), backPatchFaces);
891 updateCellSet(map(), addedCellsSet);
894 if (map().hasMotionPoints())
907 if (isType<extrudeModels::wedge>(model()))
909 changeFrontBackPatches<wedgePolyPatch>
916 else if (isType<extrudeModels::plane>(model()))
918 changeFrontBackPatches<emptyPolyPatch>
936 <<
"Cannot stitch front and back of extrusion since" 937 <<
" in 'mesh' mode (extrusion appended to mesh)." 941 Info<<
"Assuming full 360 degree axisymmetric case;" 942 <<
" stitching faces on patches " 943 << frontPatchName <<
" and " 944 << backPatchName <<
" together ..." <<
nl <<
endl;
946 if (frontPatchFaces.
size() != backPatchFaces.
size())
949 <<
"Differing number of faces on front (" 950 << frontPatchFaces.
size() <<
") and back (" 951 << backPatchFaces.
size() <<
")" 960 const word cutZoneName(
"originalCutFaceZone");
991 perfectStitcher.setRefinement
1021 updateCellSet(map(), addedCellsSet);
1024 if (map().hasMotionPoints())
1043 cellSet addedCells(mesh,
"addedCells", addedCellsSet);
1044 Info<<
"Writing added cells to cellSet " << addedCells.name()
1046 if (!addedCells.write())
1049 << addedCells.name()
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
ListType renumber(const labelUList &oldToNew, const ListType &)
Renumber the values (not the indices) of a list.
autoPtr< IOobject > clone() const
Clone.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
#define forAll(list, i)
Loop across all elements in list.
mode_t mode(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file mode.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
const word & name() const
Return name.
A class for handling file names.
errorManipArg< error, int > exit(error &err, const int errNo=1)
void reset(T *=nullptr)
If object pointer already set, delete object and set to given.
A list of keyword definitions, which are a keyword followed by any number of values (e...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
virtual const labelList & faceNeighbour() const
Return face neighbour.
const labelListList & pointEdges() const
void size(const label)
Override size to be inconsistent with allocated storage.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
static word defaultRegion
Return the default region name.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Does polyTopoChanges to remove edges. Can remove faces due to edge collapse but can not remove cells ...
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
label findPatchID(const word &patchName) const
Find patch index given a name.
label collapseEdge(triSurface &surf, const scalar minLen)
Keep collapsing all edges < minLen.
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
A bounding box defined in terms of the points at its extremities.
void transfer(HashTable< T, Key, Hash > &)
Transfer the contents of the argument table into this table.
bool insert(const Key &key)
Insert a new entry.
T & first()
Return the first element of the list.
IOdictionary systemDict(const word &dictName, const argList &args, const objectRegistry &ob, const word ®ionName=polyMesh::defaultRegion)
Initialise the NamedEnum HashTable from the static list of names.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
label size() const
Return number of elements in table.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Neighbour processor patch.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
virtual const pointField & points() const
Return raw points.
List< bool > boolList
Bool container classes.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
A list of faces which address into the list of points.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
const word & regionDir(const word ®ionName)
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
static autoPtr< extrudeModel > New(const dictionary &)
Select null constructed.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
virtual void updateMesh(const mapPolyMesh &mpm)
Update the mesh corresponding to given map.
List of mesh modifiers defining the mesh dynamics.
Adds layers of cells to outside of polyPatch. Can optionally create stand-alone extruded mesh (addToM...
A class for handling words, derived from string.
wordList names() const
Return a list of patch names.
word name() const
Return file name (part beyond last /)
virtual const labelList & faceOwner() const
Return face owner.
static const word null
An empty word.
fileName localObjectPath() const
Return complete localPath + object name.
scalar mag(const pointField &) const
Return scalar magnitude.
virtual const faceList & faces() const
Return raw faces.
void addZones(const List< pointZone *> &pz, const List< faceZone *> &fz, const List< cellZone *> &cz)
Add mesh zones.
static word controlDictName
The default control dictionary name (normally "controlDict")
bool valid() const
Return true if the autoPtr valid (ie, the pointer is set)
void addPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches.
const word & system() const
Return system name.
void removeBoundary()
Remove boundary patches.
const labelList & reverseCellMap() const
Reverse cell map.
Type gMax(const FieldField< Field, Type > &f)
const Time & time() const
Return time.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
void setInstance(const fileName &)
Set the instance for mesh files.
label size() const
Return the number of elements in the UPtrList.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
List< word > wordList
A List of words.
Calculates points shared by more than two processor patches or cyclic patches.
void setSize(const label)
Reset size of List.
label whichZone(const label objectIndex) const
Given a global object index, return the zone it is in.
const boundBox & bounds() const
Return mesh bounding box.
static bool & parRun()
Is this a parallel run?
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
string & expand(const bool allowEmpty=false)
Expand initial tildes and all occurrences of environment variables.
A collection of cell labels.
prefixOSstream Pout(cout, "Pout")
const meshFaceZones & faceZones() const
Return face zones.
label index() const
Return the index of this patch in the boundaryMesh.
Direct mesh changes based on v1.3 polyTopoChange syntax.
label start() const
Return start label of this patch in the polyMesh face list.
vector span() const
The bounding box span (from minimum to maximum)
dimensioned< scalar > mag(const dimensioned< Type > &)
static labelListList globalEdgeFaces(const polyMesh &, const globalIndex &globalFaces, const indirectPrimitivePatch &pp)
Per patch edge the pp faces (in global indices) using it. Uses.
Enum read(Istream &) const
Read a word from Istream and return the corresponding.
fileName path() const
Return directory path name (part before last /)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Mesh consisting of general polyhedral cells.
A subset of mesh faces organised as a primitive patch.
Hack of attachDetach to couple patches when they perfectly align. Does not decouple. Used by stitchMesh app. Does geometric matching.
virtual bool write(const bool write=true) const
Write using setting from DB.
A patch is a list of labels that address the faces in the global face list.
T & last()
Return the last element of the list.
Foam::argList args(argc, argv)
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
static void calcSidePatch(const polyMesh &, const globalIndex &globalFaces, const labelListList &globalEdgeFaces, const indirectPrimitivePatch &pp, labelList &sidePatchID, label &nPatches, Map< label > &nbrProcToPatch, Map< label > &patchToNbrProc)
Boundary edges get extruded into boundary faces. Determine patch.
A List with indirect addressing.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
scalar minDim() const
Smallest length/height/width dimension.
const labelList & reverseFaceMap() const
Reverse face map.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.