Go to the documentation of this file.
44 namespace functionObjects
54 void Foam::functionObjects::interfaceHeight::writePositions()
60 if (
mag(direction_) > 0.0)
62 gHat = direction_/
mag(direction_);
66 gHat =
g.value()/
mag(
g.value());
72 autoPtr<interpolation<scalar>>
87 const midPointAndFaceSet set
98 scalar hLB = set.size() ? - gHat & (locations_[li] - set[0]) : - GREAT;
99 reduce(hLB, maxOp<scalar>());
104 scalar sumLength = 0, sumLengthAlpha = 0;
105 for(label si = 0; si < set.size() - 1; ++ si)
107 if (set.segments()[si] != set.segments()[si+1])
112 const vector&
p0 = set[si], p1 = set[si+1];
113 const label c0 = set.cells()[si],
c1 = set.cells()[si+1];
114 const label f0 = set.faces()[si], f1 = set.faces()[si+1];
115 const scalar
a0 = interpolator->interpolate(
p0, c0, f0);
116 const scalar a1 = interpolator->interpolate(p1,
c1, f1);
118 const scalar l = - gHat & (p1 -
p0);
120 sumLengthAlpha += l*(
a0 + a1)/2;
123 reduce(sumLength, sumOp<scalar>());
124 reduce(sumLengthAlpha, sumOp<scalar>());
131 liquid_ ? sumLengthAlpha : sumLength - sumLengthAlpha;
132 const scalar hIL = hIB - hLB;
135 const point p = locations_[li] - gHat*hIL;
139 files(fileID::heightFile) << w << hIB << w << hIL;
140 files(fileID::positionFile) <<
'(' << w <<
p.x() << w <<
p.y()
147 files(fileID::heightFile).endl();
148 files(fileID::positionFile).endl();
169 case fileID::heightFile:
175 "Interface height above the boundary"
181 "Interface height above the location"
185 case fileID::positionFile:
187 writeHeaderValue(files(fid),
"p",
"Interface position");
194 writeCommented(files(fid),
"Location");
199 case fileID::heightFile:
200 files(fid) << w << li << w <<
' ';
202 case fileID::positionFile:
203 files(fid) << w << li << w <<
' ' << w <<
' ' <<
" ";
209 writeCommented(files(fid),
"Time");
214 case fileID::heightFile:
215 files(fid) << w <<
"hB" << w <<
"hL";
217 case fileID::positionFile:
218 files(fid) << w <<
"p" << w <<
' ' << w <<
' ' <<
" ";
239 interpolationScheme_(
"cellPoint"),
244 resetNames({
"height",
"position"});
246 writeFileHeader(fileID::heightFile);
247 writeFileHeader(fileID::positionFile);
255 dict.readIfPresent(
"alpha", alphaName_);
256 dict.readIfPresent(
"liquid", liquid_);
257 dict.readEntry(
"locations", locations_);
258 dict.readIfPresent(
"interpolationScheme", interpolationScheme_);
259 dict.readIfPresent(
"direction", direction_);
scalar mag() const
The magnitude of the bounding box span.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
virtual bool write()
Write.
virtual bool read(const dictionary &)
Read.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual bool execute()
Execute.
bool read(const char *buf, int32_t &val)
Same as readInt32.
PtrList< OFstream > & files()
Return access to the files.
static word timeName(const scalar t, const int precision=precision_)
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
#define forAll(list, i)
Loop across all elements in list.
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
static autoPtr< interpolation< scalar > > New(const word &interpolationType, const GeometricField< scalar, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
Omanip< int > valueWidth(const label offset=0) const
Return the value width when writing to stream with optional offset.
virtual bool end()
Execute at the final time-loop.
word name(const complex &c)
Return string representation of complex.
UniformDimensionedField< vector > uniformDimensionedVectorField
Istream and Ostream manipulators taking arguments.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const Type & lookupObject(const word &name, const bool recursive=false) const
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.
const uniformDimensionedVectorField & g
Vector< scalar > vector
A scalar version of the templated Vector.
functionObject base class for creating, maintaining and writing log files e.g. integrated or averaged...
static bool master(const label communicator=0)
Am I the master process.
const boundBox & bounds() const
Return mesh bounding box.
const dimensionedScalar a0
Bohr radius: default SI units: [m].
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
An Ostream manipulator taking arguments.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
const fvMesh & mesh_
Reference to the fvMesh.
defineTypeNameAndDebug(ObukhovLength, 0)
const Time & time() const
Return the top-level database.
static const Vector< scalar > zero
const volScalarField & p0
vector point
Point is a vector.
interfaceHeight(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
virtual bool write()
Write function.
virtual void writeFileHeader(const fileID fid)
Output file header information.