42 namespace functionObjects
57 Log <<
" functionObjects::" <<
type() <<
" " <<
name()
58 <<
" calculating stream-function" <<
endl;
85 label nVisitedOld = 0;
93 unitAreas /=
mag(unitAreas);
111 if (!isType<emptyPolyPatch>(
patches[patchi]))
115 if (
magSqr(
phi.boundaryField()[patchi][facei]) < SMALL)
117 const labelList& zeroPoints = bouFaces[facei];
122 forAll(zeroPoints, pointi)
124 if (visitedPoint[zeroPoints[pointi]] == 1)
133 Log <<
" Zero face: patch: " << patchi
134 <<
" face: " << facei <<
endl;
136 forAll(zeroPoints, pointi)
139 visitedPoint[zeroPoints[pointi]] = 1;
154 Log <<
" Zero flux boundary face not found. " 155 <<
"Using cell as a reference." 166 forAll(zeroPoints, pointi)
168 if (visitedPoint[zeroPoints[pointi]] == 1)
177 forAll(zeroPoints, pointi)
180 visitedPoint[zeroPoints[pointi]] = 1;
189 <<
"Cannot find initialisation face or a cell." 204 for (label facei = nInternalFaces; facei<faces.size(); facei++)
206 const labelList& curBPoints = faces[facei];
207 bool bPointFound =
false;
209 scalar currentBStream = 0.0;
210 vector currentBStreamPoint(0, 0, 0);
212 forAll(curBPoints, pointi)
215 if (visitedPoint[curBPoints[pointi]] == 1)
219 currentBStreamPoint =
points[curBPoints[pointi]];
230 forAll(curBPoints, pointi)
233 if (visitedPoint[curBPoints[pointi]] == 0)
240 !isType<emptyPolyPatch>(
patches[patchNo])
241 && !isType<symmetryPlanePolyPatch>
243 && !isType<symmetryPolyPatch>(
patches[patchNo])
244 && !isType<wedgePolyPatch>(
patches[patchNo])
252 points[curBPoints[pointi]]
253 - currentBStreamPoint;
254 edgeHat.replace(slabDir, 0);
257 vector nHat = unitAreas[facei];
259 if (edgeHat.y() > VSMALL)
261 visitedPoint[curBPoints[pointi]] = 1;
266 +
phi.boundaryField()[patchNo][faceNo]
269 else if (edgeHat.y() < -VSMALL)
271 visitedPoint[curBPoints[pointi]] = 1;
276 -
phi.boundaryField()[patchNo][faceNo]
281 if (edgeHat.x() > VSMALL)
283 visitedPoint[curBPoints[pointi]] = 1;
288 +
phi.boundaryField()[patchNo][faceNo]
291 else if (edgeHat.x() < -VSMALL)
293 visitedPoint[curBPoints[pointi]] = 1;
298 -
phi.boundaryField()[patchNo][faceNo]
312 for (label facei=0; facei<nInternalFaces; facei++)
315 const labelList& curPoints = faces[facei];
317 bool pointFound =
false;
319 scalar currentStream = 0.0;
320 point currentStreamPoint(0, 0, 0);
325 if (visitedPoint[curPoints[pointi]] == 1)
329 currentStreamPoint =
points[curPoints[pointi]];
342 if (visitedPoint[curPoints[pointi]] == 0)
345 points[curPoints[pointi]] - currentStreamPoint;
347 edgeHat.replace(slabDir, 0);
350 vector nHat = unitAreas[facei];
352 if (edgeHat.y() > VSMALL)
354 visitedPoint[curPoints[pointi]] = 1;
361 else if (edgeHat.y() < -VSMALL)
363 visitedPoint[curPoints[pointi]] = 1;
379 if (nVisited == nVisitedOld)
383 Log <<
" Exhausted a seed, looking for new seed " 384 <<
"(this is correct for multiply connected domains).";
390 nVisitedOld = nVisited;
402 return tstreamFunction;
406 bool Foam::functionObjects::streamFunction::calc()
408 const auto* phiPtr = findObject<surfaceScalarField>(fieldName_);
414 return store(resultName_, calc(
phi));
439 <<
"Case is not 2D, stream-function cannot be computed"
dimensionedScalar sign(const dimensionedScalar &ds)
defineTypeNameAndDebug(ObukhovLength, 0)
streamFunction(const word &name, const Time &runTime, const dictionary &dict)
Construct for given objectRegistry and dictionary.
List< cell > cellList
List of cell.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label nPoints() const noexcept
Number of mesh points.
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject. Registered with typeName.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const cellList & cells() const
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
Macros for easy insertion into run-time selection tables.
const word & name() const noexcept
Return the name of this functionObject.
virtual const pointField & points() const
Return raw points.
#define forAll(list, i)
Loop across all elements in list.
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
void setResultName(const word &typeName, const word &defaultArg)
Set the name of result field.
PtrList< polyPatch > polyPatchList
Store lists of polyPatch as a PtrList.
vectorField pointField
pointField is a vectorField.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
UList< face > faceUList
UList of faces.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
virtual const word & type() const =0
Runtime type information.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
label nInternalFaces() const noexcept
Number of internal faces.
virtual const faceList & faces() const
Return raw faces.
label whichPatch(const label meshFacei) const
Return patch index for a given mesh face index. Uses binary search.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
A PrimitivePatch with a SubList addressing for the faces, const reference for the point field...
vector point
Point is a vector.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Intermediate class for handling field expression function objects (e.g. blendingFactor etc...
const vectorField & faceAreas() const
const dimensionedScalar c
Speed of light in a vacuum.
const polyBoundaryMesh & patches
vector span() const
The bounding box span (from minimum to maximum)
static const Vector< label > one
const boundBox & bounds() const noexcept
Return mesh bounding box.
Field< vector > vectorField
Specialisation of Field<T> for vector.
List< label > labelList
A List of labels.
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const fvMesh & mesh_
Reference to the fvMesh.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
static constexpr const zero Zero
Global zero (0)
const Time & time_
Reference to the time database.