Go to the documentation of this file.
42 namespace functionObjects
63 Log <<
" functionObjects::" <<
type() <<
" " <<
name()
64 <<
" calculating stream-function" <<
endl;
93 label nVisitedOld = 0;
101 unitAreas /=
mag(unitAreas);
105 bool finished =
true;
119 if (!isType<emptyPolyPatch>(
patches[patchi]))
123 if (
magSqr(
phi.boundaryField()[patchi][facei]) < SMALL)
125 const labelList& zeroPoints = bouFaces[facei];
130 forAll(zeroPoints, pointi)
132 if (visitedPoint[zeroPoints[pointi]] == 1)
141 Log <<
" Zero face: patch: " << patchi
142 <<
" face: " << facei <<
endl;
144 forAll(zeroPoints, pointi)
147 visitedPoint[zeroPoints[pointi]] = 1;
162 Log <<
" Zero flux boundary face not found. "
163 <<
"Using cell as a reference."
174 forAll(zeroPoints, pointi)
176 if (visitedPoint[zeroPoints[pointi]] == 1)
185 forAll(zeroPoints, pointi)
188 visitedPoint[zeroPoints[pointi]] = 1;
197 <<
"Cannot find initialisation face or a cell."
212 for (
label facei = nInternalFaces; facei<faces.size(); facei++)
214 const labelList& curBPoints = faces[facei];
215 bool bPointFound =
false;
217 scalar currentBStream = 0.0;
218 vector currentBStreamPoint(0, 0, 0);
220 forAll(curBPoints, pointi)
223 if (visitedPoint[curBPoints[pointi]] == 1)
227 currentBStreamPoint =
points[curBPoints[pointi]];
238 forAll(curBPoints, pointi)
241 if (visitedPoint[curBPoints[pointi]] == 0)
248 !isType<emptyPolyPatch>(
patches[patchNo])
249 && !isType<symmetryPlanePolyPatch>
251 && !isType<symmetryPolyPatch>(
patches[patchNo])
252 && !isType<wedgePolyPatch>(
patches[patchNo])
260 points[curBPoints[pointi]]
261 - currentBStreamPoint;
262 edgeHat.replace(slabDir, 0);
265 vector nHat = unitAreas[facei];
267 if (edgeHat.y() > VSMALL)
269 visitedPoint[curBPoints[pointi]] = 1;
274 +
phi.boundaryField()[patchNo][faceNo]
277 else if (edgeHat.y() < -VSMALL)
279 visitedPoint[curBPoints[pointi]] = 1;
284 -
phi.boundaryField()[patchNo][faceNo]
289 if (edgeHat.x() > VSMALL)
291 visitedPoint[curBPoints[pointi]] = 1;
296 +
phi.boundaryField()[patchNo][faceNo]
299 else if (edgeHat.x() < -VSMALL)
301 visitedPoint[curBPoints[pointi]] = 1;
306 -
phi.boundaryField()[patchNo][faceNo]
320 for (
label facei=0; facei<nInternalFaces; facei++)
323 const labelList& curPoints = faces[facei];
325 bool pointFound =
false;
327 scalar currentStream = 0.0;
328 point currentStreamPoint(0, 0, 0);
333 if (visitedPoint[curPoints[pointi]] == 1)
337 currentStreamPoint =
points[curPoints[pointi]];
350 if (visitedPoint[curPoints[pointi]] == 0)
353 points[curPoints[pointi]] - currentStreamPoint;
355 edgeHat.replace(slabDir, 0);
358 vector nHat = unitAreas[facei];
360 if (edgeHat.y() > VSMALL)
362 visitedPoint[curPoints[pointi]] = 1;
369 else if (edgeHat.y() < -VSMALL)
371 visitedPoint[curPoints[pointi]] = 1;
387 if (nVisited == nVisitedOld)
391 Log <<
" Exhausted a seed, looking for new seed "
392 <<
"(this is correct for multiply connected domains).";
398 nVisitedOld = nVisited;
407 return tstreamFunction;
411 bool Foam::functionObjects::streamFunction::calc()
413 if (foundObject<surfaceScalarField>(fieldName_))
418 return store(resultName_, calc(
phi));
436 setResultName(typeName,
"phi");
438 label nD = mesh_.nGeometricD();
443 <<
"Case is not 2D, stream-function cannot be computed"
List< label > labelList
A List of labels.
vectorField pointField
pointField is a vectorField.
virtual const pointField & points() const
Return raw points.
static const Vector< label > one
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero.
label nInternalFaces() const
Number of internal faces.
const Time & time_
Reference to the time database.
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
static word timeName(const scalar t, const int precision=precision_)
const cellList & cells() const
streamFunction(const word &name, const Time &runTime, const dictionary &dict)
Construct for given objectRegistry and dictionary.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar sign(const dimensionedScalar &ds)
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
addToRunTimeSelectionTable(functionObject, add, dictionary)
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
vector span() const
The bounding box span (from minimum to maximum)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
word name(const complex &c)
Return string representation of complex.
PtrList< polyPatch > polyPatchList
container classes for polyPatch
List< cell > cellList
A List of cells.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
Vector< scalar > vector
A scalar version of the templated Vector.
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
errorManipArg< error, int > exit(error &err, const int errNo=1)
const boundBox & bounds() const
Return mesh bounding box.
defineTypeNameAndDebug(fvMeshFunctionObject, 0)
virtual const faceList & faces() const
Return raw faces.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
UList< face > faceUList
A UList of faces.
Base class for field expression function objects.
const word & name() const
Return the name of this functionObject.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual const word & type() const =0
Runtime type information.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const fvMesh & mesh_
Reference to the fvMesh.
label nPoints() const
Number of mesh points.
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
const polyBoundaryMesh & patches
const dimensionedScalar c
Speed of light in a vacuum.
PrimitivePatch< face, SubList, const pointField & > primitivePatch
Addressing for a faceList slice.
virtual ~streamFunction()
Destructor.
vector point
Point is a vector.
#define Log
Report write to Foam::Info if the local log switch is true.
const vectorField & faceAreas() const