Go to the documentation of this file.
55 void Foam::polyMesh::calcDirections()
const
65 label nEmptyPatches = 0;
66 label nWedgePatches = 0;
82 const wedgePolyPatch& wpp = refCast<const wedgePolyPatch>
88 wedgeDirVec +=
cmptMag(wpp.centreNormal());
93 reduce(nEmptyPatches, maxOp<label>());
94 reduce(nWedgePatches, maxOp<label>());
98 reduce(emptyDirVec, sumOp<vector>());
100 emptyDirVec.normalise();
104 if (emptyDirVec[cmpt] > 1
e-6)
106 solutionD_[cmpt] = -1;
110 solutionD_[cmpt] = 1;
118 geometricD_ = solutionD_;
122 reduce(wedgeDirVec, sumOp<vector>());
124 wedgeDirVec.normalise();
128 if (wedgeDirVec[cmpt] > 1
e-6)
130 geometricD_[cmpt] = -1;
134 geometricD_[cmpt] = 1;
173 time().findInstance(meshDir(),
"points"),
185 time().findInstance(meshDir(),
"faces"),
216 clearedPrimitives_(false),
240 tetBasePtIsPtr_(readTetBasePtIs()),
241 cellTreePtr_(nullptr),
299 globalMeshDataPtr_(nullptr),
301 topoChanging_(false),
303 oldPointsPtr_(nullptr)
335 boundary_.calcGeometry();
341 <<
"mesh missing boundary on one or more domains" <<
endl;
346 <<
"no points in mesh" <<
endl;
351 <<
"no cells in mesh" <<
endl;
360 Foam::polyMesh::polyMesh
424 clearedPrimitives_(
false),
439 bounds_(points_, syncPar),
443 tetBasePtIsPtr_(readTetBasePtIs()),
444 cellTreePtr_(
nullptr),
487 globalMeshDataPtr_(
nullptr),
489 topoChanging_(
false),
491 oldPointsPtr_(
nullptr)
496 const face& curFace = faces_[facei];
498 if (
min(curFace) < 0 ||
max(curFace) > points_.size())
501 <<
"Face " << facei <<
"contains vertex labels out of range: "
502 << curFace <<
" Max point index = " << points_.size()
512 Foam::polyMesh::polyMesh
575 clearedPrimitives_(
false),
590 bounds_(points_, syncPar),
594 tetBasePtIsPtr_(readTetBasePtIs()),
595 cellTreePtr_(
nullptr),
638 globalMeshDataPtr_(
nullptr),
640 topoChanging_(
false),
642 oldPointsPtr_(
nullptr)
647 const face& curFace = faces_[facei];
649 if (
min(curFace) < 0 ||
max(curFace) > points_.size())
652 <<
"Face " << facei <<
"contains vertex labels out of range: "
653 << curFace <<
" Max point index = " << points_.size()
664 const cell& curCell = cLst[celli];
666 if (
min(curCell) < 0 ||
max(curCell) > faces_.size())
669 <<
"Cell " << celli <<
"contains face labels out of range: "
670 << curCell <<
" Max face index = " << faces_.size()
680 Foam::polyMesh::polyMesh(
const IOobject& io,
const zero,
const bool syncPar)
694 const bool validBoundary
698 clearAddressing(
true);
704 points_.transfer(
points());
705 bounds_ =
boundBox(points_, validBoundary);
710 faces_.transfer(faces());
715 owner_.transfer(owner());
718 if (neighbour.valid())
720 neighbour_.transfer(neighbour());
744 const face& curFace = faces_[facei];
746 if (
min(curFace) < 0 ||
max(curFace) > points_.size())
749 <<
"Face " << facei <<
" contains vertex labels out of range: "
750 << curFace <<
" Max point index = " << points_.size()
768 boundary_.updateMesh();
771 boundary_.calcGeometry();
781 <<
"no points or no cells in mesh" <<
endl;
802 return parent().dbDir();
811 return dbDir()/meshSubDir;
817 return points_.instance();
823 return faces_.instance();
829 if (geometricD_.x() == 0)
846 if (solutionD_.x() == 0)
863 if (tetBasePtIsPtr_.empty())
868 <<
"Forcing storage of base points."
872 tetBasePtIsPtr_.reset
890 return *tetBasePtIsPtr_;
897 if (cellTreePtr_.empty())
925 return cellTreePtr_();
932 const bool validBoundary
938 <<
"boundary already exists"
946 boundary_.transfer(plist);
952 globalMeshDataPtr_.clear();
957 boundary_.updateMesh();
960 boundary_.calcGeometry();
962 boundary_.checkDefinition();
974 if (pointZones().size() || faceZones().size() || cellZones().size())
977 <<
"point, face or cell zone already exists"
984 pointZones_.setSize(pz.size());
989 pointZones_.set(pI, pz[pI]);
998 faceZones_.setSize(fz.size());
1003 faceZones_.set(fI, fz[fI]);
1012 cellZones_.setSize(cz.size());
1017 cellZones_.set(cI, cz[cI]);
1028 const bool validBoundary
1034 addPatches(plist, validBoundary);
1040 if (clearedPrimitives_)
1043 <<
"points deallocated"
1059 io.
eventNo() = points_.eventNo()+1;
1065 if (clearedPrimitives_)
1068 <<
"faces deallocated"
1090 if (oldPointsPtr_.empty())
1097 oldPointsPtr_.reset(
new pointField(points_));
1098 curMotionTimeIndex_ = time().timeIndex();
1101 return *oldPointsPtr_;
1111 <<
"Moving points for time " << time().value()
1112 <<
" index " << time().timeIndex() <<
endl;
1114 if (newPoints.size() != points_.size())
1117 <<
"Size of newPoints " << newPoints.size()
1118 <<
" does not correspond to current mesh points size "
1127 if (curMotionTimeIndex_ != time().
timeIndex())
1131 Info<<
"tmp<scalarField> polyMesh::movePoints(const pointField&) : "
1132 <<
" Storing current points for time " << time().value()
1133 <<
" index " << time().timeIndex() <<
endl;
1137 oldPointsPtr_.clear();
1138 oldPointsPtr_.reset(
new pointField(points_));
1139 curMotionTimeIndex_ = time().timeIndex();
1142 points_ = newPoints;
1144 bool moveError =
false;
1148 if (checkMeshMotion(points_,
true))
1153 <<
"Moving the mesh with given points will "
1154 <<
"invalidate the mesh." <<
nl
1155 <<
"Mesh motion should not be executed." <<
endl;
1160 points_.instance() = time().timeName();
1161 points_.eventNo() = getEvent();
1163 if (tetBasePtIsPtr_.valid())
1166 tetBasePtIsPtr_().instance() = time().timeName();
1167 tetBasePtIsPtr_().eventNo() = getEvent();
1177 if (globalMeshDataPtr_.valid())
1179 globalMeshDataPtr_().movePoints(points_);
1185 boundary_.movePoints(points_);
1187 pointZones_.movePoints(points_);
1188 faceZones_.movePoints(points_);
1189 cellZones_.movePoints(points_);
1200 cellTreePtr_.clear();
1210 meshObject::movePoints<polyMesh>(*
this);
1211 meshObject::movePoints<pointMesh>(*
this);
1213 const_cast<Time&
>(time()).functionObjects().movePoints(*
this);
1216 if (
debug && moveError)
1229 curMotionTimeIndex_ = 0;
1230 oldPointsPtr_.clear();
1236 if (globalMeshDataPtr_.empty())
1240 Pout<<
"polyMesh::globalData() const : "
1241 <<
"Constructing parallelData from processor topology"
1248 return *globalMeshDataPtr_;
1266 fileName meshFilesPath = thisDb().time().
path()/instanceDir/meshDir();
1268 rm(meshFilesPath/
"points");
1269 rm(meshFilesPath/
"faces");
1270 rm(meshFilesPath/
"owner");
1271 rm(meshFilesPath/
"neighbour");
1272 rm(meshFilesPath/
"cells");
1273 rm(meshFilesPath/
"boundary");
1274 rm(meshFilesPath/
"pointZones");
1275 rm(meshFilesPath/
"faceZones");
1276 rm(meshFilesPath/
"cellZones");
1277 rm(meshFilesPath/
"meshModifiers");
1278 rm(meshFilesPath/
"parallelData");
1281 if (
isDir(meshFilesPath/
"sets"))
1283 rmDir(meshFilesPath/
"sets");
1290 removeFiles(instance());
1314 findTetFacePt(celli,
p, tetFacei, tetPti);
1330 tetFacei = tet.
face();
1331 tetPti = tet.
tetPt();
1350 case FACE_CENTRE_TRIS:
1358 label facei = cFaces[cFacei];
1359 const face&
f = faces_[facei];
1360 const point& fc = faceCentres()[facei];
1361 bool isOwn = (owner_[facei] == celli);
1371 nextPointi =
f.nextLabel(fp);
1375 pointi =
f.nextLabel(fp);
1398 case FACE_DIAG_TRIS:
1406 label facei = cFaces[cFacei];
1407 const face&
f = faces_[facei];
1409 for (label tetPti = 1; tetPti <
f.size() - 1; tetPti++)
1434 findTetFacePt(celli,
p, tetFacei, tetPti);
1436 return tetFacei != -1;
1454 && (decompMode == FACE_DIAG_TRIS || decompMode == CELL_TETS)
1463 (void)tetBasePtIs();
1471 if (decompMode == CELL_TETS)
1480 findCellFacePt(
p, celli, tetFacei, tetPti);
1495 (void)tetBasePtIs();
1499 label celli = findNearestCell(
p);
1502 if (pointInCell(
p, celli, decompMode))
1510 for (label celli = 0; celli < nCells(); celli++)
1512 if (pointInCell(
p, celli, decompMode))
int debug
Static debugging option.
void addPatches(PtrList< polyPatch > &plist, const bool validBoundary=true)
Add boundary patches.
static autoPtr< T > New(Args &&... args)
Construct autoPtr of T with forwarding arguments.
vectorField pointField
pointField is a vectorField.
label size() const noexcept
The number of elements in table.
virtual const pointField & points() const
Return raw points.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
treeBoundBox extend(Random &rndGen, const scalar s) const
Return slightly wider bounding box.
#define InfoInFunction
Report an information message using Foam::Info.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
void resetMotion() const
Reset motion.
cellDecomposition
Enumeration defining the decomposition of the cell for.
static Vector< Cmpt > uniform(const Cmpt &s)
Return a VectorSpace with all elements = s.
A class for handling file names.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
static word defaultRegion
Return the default region name.
static std::string path(const std::string &str)
Return directory path name (part before last /)
A class for managing temporary objects.
bool upToDate(const regIOobject &) const
Return true if up-to-date with respect to given object.
static constexpr const zero Zero
Global zero (0)
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
bool typeHeaderOk(const bool checkType=true, const bool search=true, const bool verbose=true)
Read header (uses typeFilePath to find file) and check its info.
static labelList findFaceBasePts(const polyMesh &mesh, scalar tol=minTetQuality, bool report=false)
Find a suitable base point for each face for decomposition.
virtual const fileName & dbDir() const
Override the objectRegistry dbDir for a single-region case.
Standard boundBox with extra functionality for use in octree.
const Time & time() const
Return time.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
static bool & parRun()
Is this a parallel run?
bool pointInCell(const point &p, label celli) const
Return true if the point is in the cell.
label face() const
Return the face.
const fileName & facesInstance() const
Return the current instance directory for faces.
bool rm(const fileName &file)
Remove a file (or its gz equivalent), returning true if successful.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const point & max() const
Maximum describing the bounding box.
label eventNo() const
Event number at last update.
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
void removeFiles() const
Remove all files from mesh instance()
IOList< label > labelIOList
Label container classes.
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
bool pointInCell(const point &p, label celli, const cellDecomposition=CELL_TETS) const
Test if point p is in the celli.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Mesh consisting of general polyhedral cells.
const point & min() const
Minimum describing the bounding box.
#define forAll(list, i)
Loop across all elements in list.
Point centre() const
Return centre (centroid)
Registry of regIOobjects.
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
A triangle primitive used to calculate face normals and swept volumes.
const fileName & pointsInstance() const
Return the current instance directory for points.
virtual const fileName & dbDir() const
Local directory path of this objectRegistry relative to the time.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
label nCells() const
Number of mesh cells.
Encapsulation of data needed to search in/for cells. Used to find the cell containing a point (e....
Inter-processor communications stream.
virtual bool write(const bool valid=true) const
Write using setting from DB.
writeOption writeOpt() const
The write option.
triPointRef faceTri(const polyMesh &mesh) const
Return the geometry corresponding to the tri on the face for.
messageStream Info
Information stream (uses stdout - output is on the master only)
#define DebugInFunction
Report an information message using Foam::Info.
A patch is a list of labels that address the faces in the global face list.
PtrList< polyPatch > polyPatchList
container classes for polyPatch
virtual const labelList & faceOwner() const
Return face owner.
void resetPrimitives(autoPtr< pointField > &&points, autoPtr< faceList > &&faces, autoPtr< labelList > &&owner, autoPtr< labelList > &&neighbour, const labelUList &patchSizes, const labelUList &patchStarts, const bool validBoundary=true)
Reset mesh primitive data. Assumes all patch info correct.
A List of objects of type <T> with automated input and output using a compact storage....
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
label findInside(const point &) const
Find shape containing point. Only implemented for certain.
label comm() const
Return communicator used for parallel communication.
virtual const pointField & oldPoints() const
Return old points for mesh motion.
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
void findCellFacePt(const point &p, label &celli, label &tetFacei, label &tetPti) const
Find the cell, tetFacei and tetPti for point p.
virtual bool upToDatePoints(const regIOobject &io) const
Return true if io is up-to-date with points.
static tetIndices findTet(const polyMesh &mesh, label cI, const point &pt)
Find the tet decomposition of the cell containing the given point.
Cmpt cmptSum(const SphericalTensor< Cmpt > &st)
Return the sum of components of a SphericalTensor.
errorManip< error > abort(error &err)
label tetPt() const
Return the characterising tetPtI.
Vector< scalar > vector
A scalar version of the templated Vector.
fileName meshDir() const
Return the local mesh directory (dbDir()/meshSubDir)
errorManipArg< error, int > exit(error &err, const int errNo=1)
const word & headerClassName() const
Return name of the class name read from header.
void updateMesh()
Correct polyBoundaryMesh after topology update.
label nSolutionD() const
Return the number of valid solved-for dimensions in the mesh.
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
const indexedOctree< treeDataCell > & cellTree() const
Return the cell search tree.
virtual ~polyMesh()
Destructor.
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
label findCell(const point &p, const cellDecomposition=CELL_TETS) const
Find cell enclosing this location and return index.
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
bool rmDir(const fileName &directory, const bool silent=false)
Remove a directory and its contents (optionally silencing warnings)
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
static label worldComm
Default communicator (all processors)
readOption readOpt() const
The read option.
const dimensionedScalar e
Elementary charge.
label nPoints() const
Number of mesh points.
tmp< scalarField > movePoints(const pointField &p, const pointField &oldP)
Move points, returns volumes swept by faces in motion.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
A bounding box defined in terms of min/max extrema points.
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
A face is a list of labels corresponding to mesh vertices.
const globalMeshData & globalData() const
Return parallel info.
vector areaNormal() const
The area normal - with magnitude equal to area of triangle.
void findTetFacePt(const label celli, const point &p, label &tetFacei, label &tetPti) const
Find the tetFacei and tetPti for point p in celli.
defineTypeNameAndDebug(combustionModel, 0)
#define WarningInFunction
Report a warning using Foam::Warning.
A cell is defined as a list of faces with extra functionality.
void addZones(const List< pointZone * > &pz, const List< faceZone * > &fz, const List< cellZone * > &cz)
Add mesh zones.
virtual void setUpToDatePoints(regIOobject &io) const
Set io to be up-to-date with points.
virtual const labelList & faceNeighbour() const
Return face neighbour.
static constexpr direction nComponents
Number of components in this vector space.
bool isDir(const fileName &name, const bool followLink=true)
Does the name exist as a DIRECTORY in the file system?
A class representing the concept of 0 (zero) that can be used to avoid manipulating objects known to ...
const vectorField & faceAreas() const
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Cell-face mesh analysis engine.