39 namespace functionObjects
60 Log <<
" functionObjects::" <<
type() <<
" " <<
name()
61 <<
" calculating steam-function" <<
endl;
70 scalar thickness =
vector(slabNormal) & mesh_.bounds().span();
74 tmp<pointScalarField> tstreamFunction
91 forAll(visitedPoint, pointi)
93 visitedPoint[pointi] = 0;
96 label nVisitedOld = 0;
101 label nInternalFaces = mesh_.nInternalFaces();
104 unitAreas /=
mag(unitAreas);
108 bool finished =
true;
122 if (!isType<emptyPolyPatch>(patches[patchi]))
131 const labelList& zeroPoints = bouFaces[facei];
136 forAll(zeroPoints, pointi)
138 if (visitedPoint[zeroPoints[pointi]] == 1)
147 Log <<
" Zero face: patch: " << patchi
148 <<
" face: " << facei <<
endl;
150 forAll(zeroPoints, pointi)
152 streamFunction[zeroPoints[pointi]] = 0;
153 visitedPoint[zeroPoints[pointi]] = 1;
168 Log <<
" Zero flux boundary face not found. " 169 <<
"Using cell as a reference." 176 labelList zeroPoints = c[cI].labels(mesh_.faces());
180 forAll(zeroPoints, pointi)
182 if (visitedPoint[zeroPoints[pointi]] == 1)
191 forAll(zeroPoints, pointi)
193 streamFunction[zeroPoints[pointi]] = 0.0;
194 visitedPoint[zeroPoints[pointi]] = 1;
203 <<
"Cannot find initialisation face or a cell." 218 for (
label facei = nInternalFaces; facei<faces.size(); facei++)
220 const labelList& curBPoints = faces[facei];
221 bool bPointFound =
false;
223 scalar currentBStream = 0.0;
224 vector currentBStreamPoint(0, 0, 0);
226 forAll(curBPoints, pointi)
229 if (visitedPoint[curBPoints[pointi]] == 1)
232 currentBStream = streamFunction[curBPoints[pointi]];
233 currentBStreamPoint = points[curBPoints[pointi]];
244 forAll(curBPoints, pointi)
247 if (visitedPoint[curBPoints[pointi]] == 0)
250 mesh_.boundaryMesh().whichPatch(facei);
254 !isType<emptyPolyPatch>(patches[patchNo])
255 && !isType<symmetryPlanePolyPatch>
257 && !isType<symmetryPolyPatch>(patches[patchNo])
258 && !isType<wedgePolyPatch>(patches[patchNo])
262 mesh_.boundaryMesh()[patchNo]
266 points[curBPoints[pointi]]
267 - currentBStreamPoint;
269 edgeHat /=
mag(edgeHat);
271 vector nHat = unitAreas[facei];
273 if (edgeHat.y() > vSmall)
275 visitedPoint[curBPoints[pointi]] = 1;
278 streamFunction[curBPoints[pointi]] =
280 + phi.boundaryField()[patchNo][faceNo]
283 else if (edgeHat.y() < -vSmall)
285 visitedPoint[curBPoints[pointi]] = 1;
288 streamFunction[curBPoints[pointi]] =
290 - phi.boundaryField()[patchNo][faceNo]
295 if (edgeHat.x() > vSmall)
297 visitedPoint[curBPoints[pointi]] = 1;
300 streamFunction[curBPoints[pointi]] =
302 + phi.boundaryField()[patchNo][faceNo]
305 else if (edgeHat.x() < -vSmall)
307 visitedPoint[curBPoints[pointi]] = 1;
310 streamFunction[curBPoints[pointi]] =
312 - phi.boundaryField()[patchNo][faceNo]
326 for (
label facei=0; facei<nInternalFaces; facei++)
329 const labelList& curPoints = faces[facei];
331 bool pointFound =
false;
333 scalar currentStream = 0.0;
334 point currentStreamPoint(0, 0, 0);
339 if (visitedPoint[curPoints[pointi]] == 1)
342 currentStream = streamFunction[curPoints[pointi]];
343 currentStreamPoint = points[curPoints[pointi]];
356 if (visitedPoint[curPoints[pointi]] == 0)
359 points[curPoints[pointi]] - currentStreamPoint;
362 edgeHat /=
mag(edgeHat);
364 vector nHat = unitAreas[facei];
366 if (edgeHat.y() > vSmall)
368 visitedPoint[curPoints[pointi]] = 1;
371 streamFunction[curPoints[pointi]] =
373 + phi[facei]*
sign(nHat.x());
375 else if (edgeHat.y() < -vSmall)
377 visitedPoint[curPoints[pointi]] = 1;
380 streamFunction[curPoints[pointi]] =
382 - phi[facei]*
sign(nHat.x());
393 if (nVisited == nVisitedOld)
397 Log <<
" Exhausted a seed, looking for new seed " 398 <<
"(this is correct for multiply connected domains).";
404 nVisitedOld = nVisited;
410 streamFunction /= thickness;
411 streamFunction.boundaryFieldRef() = 0.0;
413 return tstreamFunction;
417 bool Foam::functionObjects::streamFunction::calc()
419 if (foundObject<surfaceScalarField>(fieldName_))
424 return store(resultName_, calc(phi));
444 setResultName(
"streamFunction",
"phi");
446 label nD = mesh_.nGeometricD();
451 <<
"Case is not 2D, stream-function cannot be computed"
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
PtrList< polyPatch > polyPatchList
container classes for polyPatch
dimensionedScalar sign(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
streamFunction(const word &name, const Time &runTime, const dictionary &dict)
Construct for given objectRegistry and dictionary.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
virtual ~streamFunction()
Destructor.
A list of keyword definitions, which are a keyword followed by any number of values (e...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Vector< scalar > vector
A scalar version of the templated Vector.
void replace(const direction, const Cmpt &)
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Macros for easy insertion into run-time selection tables.
static const pointMesh & New(const polyMesh &mesh)
vectorField pointField
pointField is a vectorField.
A class for handling words, derived from string.
List< label > labelList
A List of labels.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
PrimitivePatch< face, SubList, const pointField & > primitivePatch
Addressing for a faceList slice.
word name(const complex &)
Return a string representation of a complex.
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Internal & ref()
Return a reference to the dimensioned internal field.
vector point
Point is a vector.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
#define Log
Report write to Foam::Info if the local log switch is true.
static const Vector< label > one
dimensioned< scalar > mag(const dimensioned< Type > &)
Field< vector > vectorField
Specialisation of Field<T> for vector.
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
A class for managing temporary objects.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
List< cell > cellList
list of cells
addToRunTimeSelectionTable(functionObject, add, dictionary)