39 namespace functionObjects
60 Log <<
" functionObjects::" <<
type() <<
" " <<
name()
61 <<
" calculating steam-function" <<
endl;
74 tmp<pointScalarField> tstreamFunction
86 forAll(visitedPoint, pointi)
88 visitedPoint[pointi] = 0;
91 label nVisitedOld = 0;
99 unitAreas /=
mag(unitAreas);
103 bool finished =
true;
126 const labelList& zeroPoints = bouFaces[facei];
131 forAll(zeroPoints, pointi)
133 if (visitedPoint[zeroPoints[pointi]] == 1)
143 <<
" face: " << facei <<
endl;
145 forAll(zeroPoints, pointi)
148 visitedPoint[zeroPoints[pointi]] = 1;
163 Log <<
" Zero flux boundary face not found. "
164 <<
"Using cell as a reference."
175 forAll(zeroPoints, pointi)
177 if (visitedPoint[zeroPoints[pointi]] == 1)
186 forAll(zeroPoints, pointi)
189 visitedPoint[zeroPoints[pointi]] = 1;
198 <<
"Cannot find initialisation face or a cell."
213 for (
label facei = nInternalFaces; facei<faces.size(); facei++)
215 const labelList& curBPoints = faces[facei];
216 bool bPointFound =
false;
218 scalar currentBStream = 0.0;
219 vector currentBStreamPoint(0, 0, 0);
221 forAll(curBPoints, pointi)
224 if (visitedPoint[curBPoints[pointi]] == 1)
228 currentBStreamPoint =
points[curBPoints[pointi]];
239 forAll(curBPoints, pointi)
242 if (visitedPoint[curBPoints[pointi]] == 0)
249 !isType<emptyPolyPatch>(
patches[patchNo])
250 && !isType<symmetryPlanePolyPatch>
252 && !isType<symmetryPolyPatch>(
patches[patchNo])
253 && !isType<wedgePolyPatch>(
patches[patchNo])
261 points[curBPoints[pointi]]
262 - currentBStreamPoint;
263 edgeHat.replace(slabDir, 0);
264 edgeHat /=
mag(edgeHat);
266 vector nHat = unitAreas[facei];
268 if (edgeHat.y() > vSmall)
270 visitedPoint[curBPoints[pointi]] = 1;
275 + phi.boundaryField()[patchNo][faceNo]
278 else if (edgeHat.y() < -vSmall)
280 visitedPoint[curBPoints[pointi]] = 1;
285 - phi.boundaryField()[patchNo][faceNo]
290 if (edgeHat.x() > vSmall)
292 visitedPoint[curBPoints[pointi]] = 1;
297 + phi.boundaryField()[patchNo][faceNo]
300 else if (edgeHat.x() < -vSmall)
302 visitedPoint[curBPoints[pointi]] = 1;
307 - phi.boundaryField()[patchNo][faceNo]
321 for (
label facei=0; facei<nInternalFaces; facei++)
324 const labelList& curPoints = faces[facei];
326 bool pointFound =
false;
328 scalar currentStream = 0.0;
329 point currentStreamPoint(0, 0, 0);
334 if (visitedPoint[curPoints[pointi]] == 1)
338 currentStreamPoint =
points[curPoints[pointi]];
351 if (visitedPoint[curPoints[pointi]] == 0)
354 points[curPoints[pointi]] - currentStreamPoint;
356 edgeHat.replace(slabDir, 0);
357 edgeHat /=
mag(edgeHat);
359 vector nHat = unitAreas[facei];
361 if (edgeHat.y() > vSmall)
363 visitedPoint[curPoints[pointi]] = 1;
368 + phi[facei]*
sign(nHat.x());
370 else if (edgeHat.y() < -vSmall)
372 visitedPoint[curPoints[pointi]] = 1;
377 - phi[facei]*
sign(nHat.x());
388 if (nVisited == nVisitedOld)
392 Log <<
" Exhausted a seed, looking for new seed "
393 <<
"(this is correct for multiply connected domains).";
399 nVisitedOld = nVisited;
408 return tstreamFunction;
412 bool Foam::functionObjects::streamFunction::calc()
414 if (foundObject<surfaceScalarField>(fieldName_))
419 return store(resultName_, calc(phi));
423 cannotFindObject<surfaceScalarField>(fieldName_);
446 <<
"Case is not 2D, stream-function cannot be computed"
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
static pointMesh & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
static tmp< GeometricField< Type, pointPatchField, pointMesh > > New(const word &name, const Internal &, const PtrList< pointPatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
vector span() const
The bounding box span (from minimum to maximum)
A list of keyword definitions, which are a keyword followed by any number of values (e....
Abstract base-class for Time/database functionObjects.
virtual const word & type() const =0
Runtime type information.
const word & name() const
Return the name of this functionObject.
const fvMesh & mesh_
Reference to the fvMesh.
This function object calculates and outputs the stream-function as a pointScalarField.
virtual ~streamFunction()
Destructor.
streamFunction(const word &name, const Time &runTime, const dictionary &dict)
Construct for given objectRegistry and dictionary.
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 polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
virtual const pointField & points() const
Return raw points.
const boundBox & bounds() const
Return mesh bounding box.
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
label nInternalFaces() const
const vectorField & faceAreas() const
const cellList & cells() const
A class for managing temporary objects.
A class for handling words, derived from string.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const fvPatchList & patches
#define Log
Report write to Foam::Info if the local log switch is true.
const dimensionedScalar c
Speed of light in a vacuum.
defineTypeNameAndDebug(adjustTimeStepToCombustion, 0)
addToRunTimeSelectionTable(functionObject, adjustTimeStepToCombustion, dictionary)
errorManipArg< error, int > exit(error &err, const int errNo=1)
List< label > labelList
A List of labels.
dimensionedScalar sign(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
List< cell > cellList
list of cells
Ostream & endl(Ostream &os)
Add newline and flush stream.
word name(const bool)
Return a word representation of a bool.
SurfaceField< scalar > surfaceScalarField
vectorField pointField
pointField is a vectorField.
vector point
Point is a vector.
PointField< scalar > pointScalarField
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
Addressing for a faceList slice.
Vector< scalar > vector
A scalar version of the templated Vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
PtrList< polyPatch > polyPatchList
container classes for polyPatch
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.