Go to the documentation of this file.
45 Foam::directions::directionTypeNames_
47 { directionType::TAN1,
"tan1" },
48 { directionType::TAN2,
"tan2" },
49 { directionType::NORMAL,
"normal" },
55 void Foam::directions::writeOBJ(Ostream& os,
const point& pt)
57 os <<
"v " << pt.x() <<
' ' << pt.y() <<
' ' << pt.z() <<
endl;
61 void Foam::directions::writeOBJ
72 os <<
"l " << vertI + 1 <<
' ' << vertI + 2 <<
endl;
78 void Foam::directions::writeOBJ
80 const fileName& fName,
81 const primitiveMesh&
mesh,
85 Pout<<
"Writing cell info to " << fName <<
" as vectors at the cellCentres"
88 OFstream xDirStream(fName);
97 scalar minDist = GREAT;
106 scalar scale = 0.5*minDist;
108 writeOBJ(xDirStream, ctr, ctr + scale*dirs[celli], vertI);
113 void Foam::directions::check2D
115 const twoDPointCorrector* correct2DPtr,
121 if (
mag(correct2DPtr->planeNormal() & vec) > 1
e-6)
124 <<
"is not normal to plane defined in dynamicMeshDict."
126 <<
"Either make case 3D or adjust vector."
135 const polyMesh&
mesh,
144 List<directionInfo> changedFacesInfo(pp.size());
150 label meshFacei = pp.start() + patchFacei;
154 if (!hexMatcher().
isA(
mesh, celli))
157 <<
"useHexTopology specified but cell " << celli
158 <<
" on face " << patchFacei <<
" of patch " << pp.name()
162 const vector& cutDir = ppField[patchFacei];
178 changedFaces[patchFacei] = meshFacei;
179 changedFacesInfo[patchFacei] =
191 changedFaces[patchFacei] = pp.start() + patchFacei;
192 changedFacesInfo[patchFacei] =
201 MeshWave<directionInfo> directionCalc
209 const List<directionInfo>& cellInfo = directionCalc.allCellInfo();
219 label index = cellInfo[celli].index();
225 <<
"Cell " << celli <<
" never visited to determine "
226 <<
"local coordinate system" <<
endl
227 <<
"Using direction " << defaultDir <<
" instead" <<
endl;
229 dirField[celli] = defaultDir;
233 else if (index == -2)
236 dirField[celli] = cellInfo[celli].n();
240 else if (index == -1)
243 <<
"Illegal index " << index <<
endl
255 reduce(nGeom, sumOp<label>());
256 reduce(nTopo, sumOp<label>());
257 reduce(nUnset, sumOp<label>());
259 Info<<
"Calculated local coords for " << defaultDir
261 <<
" Geometric cut cells : " << nGeom <<
endl
262 <<
" Topological cut cells : " << nTopo <<
endl
263 <<
" Unset cells : " << nUnset <<
endl
272 Foam::directions::directions
286 bool wantNormal =
false;
287 bool wantTan1 =
false;
288 bool wantTan2 =
false;
291 if (coordSystem !=
"fieldBased")
293 for (
const word& wantedName : wantedDirs)
297 if (wantedDir == NORMAL)
301 else if (wantedDir == TAN1)
305 else if (wantedDir == TAN2)
313 if (coordSystem ==
"global")
318 check2D(correct2DPtr, tan1);
321 check2D(correct2DPtr, tan2);
325 Info<<
"Global Coordinate system:" <<
endl
326 <<
" normal : " << normal <<
endl
327 <<
" tan1 : " << tan1 <<
endl
328 <<
" tan2 : " << tan2
344 else if (coordSystem ==
"patchLocal")
348 const word patchName(patchDict.
get<
word>(
"patch"));
355 <<
"Cannot find patch "
372 <<
"Discarding user specified tan1 since 2D case." <<
endl
373 <<
"Recalculated tan1 from face normal and planeNormal as "
377 const bool useTopo(
dict.
get<
bool>(
"useHexTopology"));
382 if (wantNormal || wantTan2)
396 this->operator[](nDirs++) = normalDirs;
400 if (wantTan1 || wantTan2)
415 this->operator[](nDirs++) = tan1Dirs;
422 this->operator[](nDirs++) = tan2Dirs;
425 else if (coordSystem ==
"fieldBased")
429 operator[](nDirs++) =
445 <<
"Unknown coordinate system "
446 << coordSystem <<
endl
447 <<
"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.
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)
Generic templated field type.
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.
void resize(const label newSize)
Adjust allocated size of list.
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.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const dimensionedScalar e
Elementary charge.
const Field< point_type > & faceNormals() const
Return face unit normals for patch.
const globalMeshData & globalData() const
Return parallel info.
vector point
Point is a vector.
#define WarningInFunction
Report a warning using Foam::Warning.