42namespace functionObjects
54 const surfaceScalarField&
phi
57 Log <<
" functionObjects::" <<
type() <<
" " <<
name()
58 <<
" calculating stream-function" <<
endl;
64 & Vector<label>(Vector<label>::X, Vector<label>::Y, Vector<label>::Z)
87 label nVisitedOld = 0;
95 unitAreas /= mag(unitAreas);
113 if (!isType<emptyPolyPatch>(
patches[patchi]))
119 const labelList& zeroPoints = bouFaces[facei];
124 forAll(zeroPoints, pointi)
126 if (visitedPoint[zeroPoints[pointi]] == 1)
135 Log <<
" Zero face: patch: " << patchi
136 <<
" face: " << facei <<
endl;
138 forAll(zeroPoints, pointi)
140 streamFunction[zeroPoints[pointi]] = 0;
141 visitedPoint[zeroPoints[pointi]] = 1;
156 Log <<
" Zero flux boundary face not found. "
157 <<
"Using cell as a reference."
168 forAll(zeroPoints, pointi)
170 if (visitedPoint[zeroPoints[pointi]] == 1)
179 forAll(zeroPoints, pointi)
181 streamFunction[zeroPoints[pointi]] = 0.0;
182 visitedPoint[zeroPoints[pointi]] = 1;
191 <<
"Cannot find initialisation face or a cell."
206 for (label facei = nInternalFaces; facei<faces.size(); facei++)
208 const labelList& curBPoints = faces[facei];
209 bool bPointFound =
false;
211 scalar currentBStream = 0.0;
212 vector currentBStreamPoint(0, 0, 0);
214 forAll(curBPoints, pointi)
217 if (visitedPoint[curBPoints[pointi]] == 1)
220 currentBStream = streamFunction[curBPoints[pointi]];
221 currentBStreamPoint =
points[curBPoints[pointi]];
232 forAll(curBPoints, pointi)
235 if (visitedPoint[curBPoints[pointi]] == 0)
242 !isType<emptyPolyPatch>(
patches[patchNo])
243 && !isType<symmetryPlanePolyPatch>
245 && !isType<symmetryPolyPatch>(
patches[patchNo])
246 && !isType<wedgePolyPatch>(
patches[patchNo])
254 points[curBPoints[pointi]]
255 - currentBStreamPoint;
256 edgeHat.replace(slabDir, 0);
259 vector nHat = unitAreas[facei];
261 if (edgeHat.y() > VSMALL)
263 visitedPoint[curBPoints[pointi]] = 1;
266 streamFunction[curBPoints[pointi]] =
271 else if (edgeHat.y() < -VSMALL)
273 visitedPoint[curBPoints[pointi]] = 1;
276 streamFunction[curBPoints[pointi]] =
283 if (edgeHat.x() > VSMALL)
285 visitedPoint[curBPoints[pointi]] = 1;
288 streamFunction[curBPoints[pointi]] =
293 else if (edgeHat.x() < -VSMALL)
295 visitedPoint[curBPoints[pointi]] = 1;
298 streamFunction[curBPoints[pointi]] =
314 for (label facei=0; facei<nInternalFaces; facei++)
317 const labelList& curPoints = faces[facei];
319 bool pointFound =
false;
321 scalar currentStream = 0.0;
322 point currentStreamPoint(0, 0, 0);
327 if (visitedPoint[curPoints[pointi]] == 1)
330 currentStream = streamFunction[curPoints[pointi]];
331 currentStreamPoint =
points[curPoints[pointi]];
344 if (visitedPoint[curPoints[pointi]] == 0)
347 points[curPoints[pointi]] - currentStreamPoint;
349 edgeHat.replace(slabDir, 0);
352 vector nHat = unitAreas[facei];
354 if (edgeHat.y() > VSMALL)
356 visitedPoint[curPoints[pointi]] = 1;
359 streamFunction[curPoints[pointi]] =
363 else if (edgeHat.y() < -VSMALL)
365 visitedPoint[curPoints[pointi]] = 1;
368 streamFunction[curPoints[pointi]] =
381 if (nVisited == nVisitedOld)
385 Log <<
" Exhausted a seed, looking for new seed "
386 <<
"(this is correct for multiply connected domains).";
392 nVisitedOld = nVisited;
398 streamFunction /= thickness;
399 streamFunction.boundaryFieldRef() = 0.0;
401 return tstreamFunction;
405bool Foam::functionObjects::streamFunction::calc()
407 if (foundObject<surfaceScalarField>(fieldName_))
412 return store(resultName_, calc(
phi));
437 <<
"Case is not 2D, stream-function cannot be computed"
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const dimensionSet & dimensions() const
Return dimensions.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
vector span() const
The bounding box span (from minimum to maximum)
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Abstract base-class for Time/database function objects.
const word & name() const noexcept
Return the name of this functionObject.
virtual const word & type() const =0
Runtime type information.
Intermediate class for handling field expression function objects (e.g. blendingFactor etc....
void setResultName(const word &typeName, const word &defaultArg)
Set the name of result field.
const fvMesh & mesh_
Reference to the fvMesh.
Computes the stream function (i.e. https://w.wiki/Ncm).
const Time & time_
Reference to the time database.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
virtual const faceList & faces() const
Return raw faces.
const boundBox & bounds() const
Return mesh bounding box.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
label nInternalFaces() const noexcept
Number of internal faces.
label nPoints() const noexcept
Number of mesh points.
const vectorField & faceAreas() const
const cellList & cells() const
A class for managing temporary objects.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
const polyBoundaryMesh & patches
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const dimensionedScalar c
Speed of light in a vacuum.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< polyPatch > polyPatchList
Store lists of polyPatch as a PtrList.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
dimensionedScalar sign(const dimensionedScalar &ds)
List< label > labelList
A List of labels.
List< cell > cellList
A List of cells.
vectorField pointField
pointField is a vectorField.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
vector point
Point is a vector.
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
UList< face > faceUList
A UList of faces.
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define forAll(list, i)
Loop across all elements in list.