Go to the documentation of this file.
44 Foam::directions::directionTypeNames_
46 { directionType::TAN1,
"tan1" },
47 { directionType::TAN2,
"tan2" },
48 { directionType::NORMAL,
"normal" },
54 void Foam::directions::writeOBJ(Ostream& os,
const point& pt)
56 os <<
"v " << pt.x() <<
' ' << pt.y() <<
' ' << pt.z() <<
endl;
60 void Foam::directions::writeOBJ
71 os <<
"l " << vertI + 1 <<
' ' << vertI + 2 <<
endl;
77 void Foam::directions::writeOBJ
79 const fileName& fName,
80 const primitiveMesh&
mesh,
84 Pout<<
"Writing cell info to " << fName <<
" as vectors at the cellCentres"
87 OFstream xDirStream(fName);
96 scalar minDist = GREAT;
105 scalar scale = 0.5*minDist;
107 writeOBJ(xDirStream, ctr, ctr + scale*dirs[celli], vertI);
112 void Foam::directions::check2D
114 const twoDPointCorrector* correct2DPtr,
120 if (
mag(correct2DPtr->planeNormal() & vec) > 1
e-6)
123 <<
"is not normal to plane defined in dynamicMeshDict."
125 <<
"Either make case 3D or adjust vector."
134 const polyMesh&
mesh,
143 List<directionInfo> changedFacesInfo(pp.size());
149 label meshFacei = pp.start() + patchFacei;
153 if (!hexMatcher().
isA(
mesh, celli))
156 <<
"useHexTopology specified but cell " << celli
157 <<
" on face " << patchFacei <<
" of patch " << pp.name()
161 const vector& cutDir = ppField[patchFacei];
177 changedFaces[patchFacei] = meshFacei;
178 changedFacesInfo[patchFacei] =
190 changedFaces[patchFacei] = pp.start() + patchFacei;
191 changedFacesInfo[patchFacei] =
200 MeshWave<directionInfo> directionCalc
208 const List<directionInfo>& cellInfo = directionCalc.allCellInfo();
218 label index = cellInfo[celli].index();
224 <<
"Cell " << celli <<
" never visited to determine "
225 <<
"local coordinate system" <<
endl
226 <<
"Using direction " << defaultDir <<
" instead" <<
endl;
228 dirField[celli] = defaultDir;
232 else if (index == -2)
235 dirField[celli] = cellInfo[celli].n();
239 else if (index == -1)
242 <<
"Illegal index " << index <<
endl
254 reduce(nGeom, sumOp<label>());
255 reduce(nTopo, sumOp<label>());
256 reduce(nUnset, sumOp<label>());
258 Info<<
"Calculated local coords for " << defaultDir
260 <<
" Geometric cut cells : " << nGeom <<
endl
261 <<
" Topological cut cells : " << nTopo <<
endl
262 <<
" Unset cells : " << nUnset <<
endl
271 Foam::directions::directions
283 bool wantNormal =
false;
284 bool wantTan1 =
false;
285 bool wantTan2 =
false;
288 if (coordSystem !=
"fieldBased")
290 for (
const word& wantedName : wantedDirs)
294 if (wantedDir == NORMAL)
298 else if (wantedDir == TAN1)
302 else if (wantedDir == TAN2)
310 if (coordSystem ==
"global")
315 check2D(correct2DPtr, tan1);
318 check2D(correct2DPtr, tan2);
322 Info<<
"Global Coordinate system:" <<
endl
323 <<
" normal : " << normal <<
endl
324 <<
" tan1 : " << tan1 <<
endl
325 <<
" tan2 : " << tan2
341 else if (coordSystem ==
"patchLocal")
345 const word patchName(patchDict.
get<
word>(
"patch"));
352 <<
"Cannot find patch "
369 <<
"Discarding user specified tan1 since 2D case." <<
endl
370 <<
"Recalculated tan1 from face normal and planeNormal as "
374 const bool useTopo(
dict.
get<
bool>(
"useHexTopology"));
379 if (wantNormal || wantTan2)
393 this->operator[](nDirs++) = normalDirs;
397 if (wantTan1 || wantTan2)
412 this->operator[](nDirs++) = tan1Dirs;
419 this->operator[](nDirs++) = tan2Dirs;
422 else if (coordSystem ==
"fieldBased")
426 operator[](nDirs++) =
442 <<
"Unknown coordinate system "
443 << coordSystem <<
endl
444 <<
"Known types are global, patchLocal and fieldBased"
List< label > labelList
A List of labels.
IOField< vector > vectorIOField
vectorField with IO.
static label edgeToFaceIndex(const primitiveMesh &mesh, const label celli, const label facei, const label edgeI)
Given edge on hex cell find corresponding edge on face. Is either.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
A class for handling words, derived from Foam::string.
const vector & planeNormal() const
Return plane normal.
A class for managing temporary objects.
const fileName & instance() const
label nTotalCells() const
Return total number of cells in decomposed mesh.
directionType
Enumeration listing the possible coordinate directions.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
prefixOSstream Pout
An Ostream wrapper for parallel output to std::cout.
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.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< word > wordList
A List of words.
Class applies a two-dimensional correction to mesh motion point field.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
messageStream Info
Information stream (uses stdout - output is on the master only)
A patch is a list of labels that address the faces in the global face list.
virtual const labelList & faceOwner() const
Return face owner.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
const labelListList & cellCells() const
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
ITstream & lookup(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
errorManip< error > abort(error &err)
Vector< scalar > vector
A scalar version of the templated Vector.
errorManipArg< error, int > exit(error &err, const int errNo=1)
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
const TargetType * isA(const Type &t)
Check if dynamic_cast to TargetType is possible.
const vectorField & cellCentres() const
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const Field< PointType > & faceNormals() const
Return face unit normals for patch.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
const globalMeshData & globalData() const
Return parallel info.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.