interfaceHeight.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2017-2019 OpenFOAM Foundation
9  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "interfaceHeight.H"
30 #include "fvMesh.H"
31 #include "interpolation.H"
32 #include "IOmanip.H"
33 #include "meshSearch.H"
34 #include "midPointAndFaceSet.H"
35 #include "Time.H"
37 #include "volFields.H"
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 namespace Foam
43 {
44 namespace functionObjects
45 {
46  defineTypeNameAndDebug(interfaceHeight, 0);
47  addToRunTimeSelectionTable(functionObject, interfaceHeight, dictionary);
48 }
49 }
50 
51 
52 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 
54 void Foam::functionObjects::interfaceHeight::writePositions()
55 {
58  vector gHat = vector::zero;
59 
60  if (mag(direction_) > 0.0)
61  {
62  gHat = direction_/mag(direction_);
63  }
64  else
65  {
66  gHat = g.value()/mag(g.value());
67  }
68 
69  const volScalarField& alpha =
70  mesh_.lookupObject<volScalarField>(alphaName_);
71 
72  autoPtr<interpolation<scalar>>
73  interpolator
74  (
75  interpolation<scalar>::New(interpolationScheme_, alpha)
76  );
77 
78  if (Pstream::master())
79  {
80  files(fileID::heightFile) << mesh_.time().timeName() << tab;
81  files(fileID::positionFile) << mesh_.time().timeName() << tab;
82  }
83 
84  forAll(locations_, li)
85  {
86  // Create a set along a ray projected in the direction of gravity
87  const midPointAndFaceSet set
88  (
89  "",
90  mesh_,
91  meshSearch(mesh_),
92  "xyz",
93  locations_[li] + gHat*mesh_.bounds().mag(),
94  locations_[li] - gHat*mesh_.bounds().mag()
95  );
96 
97  // Find the height of the location above the boundary
98  scalar hLB = set.size() ? - gHat & (locations_[li] - set[0]) : - GREAT;
99  reduce(hLB, maxOp<scalar>());
100 
101  // Calculate the integrals of length and length*alpha along the sampling
102  // line. The latter is equal to the equivalent length with alpha equal
103  // to one.
104  scalar sumLength = 0, sumLengthAlpha = 0;
105  for(label si = 0; si < set.size() - 1; ++ si)
106  {
107  if (set.segments()[si] != set.segments()[si+1])
108  {
109  continue;
110  }
111 
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);
117 
118  const scalar l = - gHat & (p1 - p0);
119  sumLength += l;
120  sumLengthAlpha += l*(a0 + a1)/2;
121  }
122 
123  reduce(sumLength, sumOp<scalar>());
124  reduce(sumLengthAlpha, sumOp<scalar>());
125 
126  // Write out
127  if (Pstream::master())
128  {
129  // Interface heights above the boundary and location
130  const scalar hIB =
131  liquid_ ? sumLengthAlpha : sumLength - sumLengthAlpha;
132  const scalar hIL = hIB - hLB;
133 
134  // Position of the interface
135  const point p = locations_[li] - gHat*hIL;
136 
137  const Foam::Omanip<int> w = valueWidth(1);
138 
139  files(fileID::heightFile) << w << hIB << w << hIL;
140  files(fileID::positionFile) << '(' << w << p.x() << w << p.y()
141  << valueWidth() << p.z() << ") ";
142  }
143  }
144 
145  if (Pstream::master())
146  {
147  files(fileID::heightFile).endl();
148  files(fileID::positionFile).endl();
149  }
150 }
151 
152 
153 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
154 
156 {
157  forAll(locations_, li)
158  {
159  writeHeaderValue
160  (
161  files(fid),
162  "Location " + Foam::name(li),
163  locations_[li]
164  );
165  }
166 
167  switch (fileID(fid))
168  {
169  case fileID::heightFile:
170  {
171  writeHeaderValue
172  (
173  files(fid),
174  "hB",
175  "Interface height above the boundary"
176  );
177  writeHeaderValue
178  (
179  files(fid),
180  "hL",
181  "Interface height above the location"
182  );
183  break;
184  }
185  case fileID::positionFile:
186  {
187  writeHeaderValue(files(fid), "p", "Interface position");
188  break;
189  }
190  }
191 
192  const Foam::Omanip<int> w = valueWidth(1);
193 
194  writeCommented(files(fid), "Location");
195  forAll(locations_, li)
196  {
197  switch (fid)
198  {
199  case fileID::heightFile:
200  files(fid) << w << li << w << ' ';
201  break;
202  case fileID::positionFile:
203  files(fid) << w << li << w << ' ' << w << ' ' << " ";
204  break;
205  }
206  }
207  files(fid).endl();
208 
209  writeCommented(files(fid), "Time");
210  forAll(locations_, li)
211  {
212  switch (fid)
213  {
214  case fileID::heightFile:
215  files(fid) << w << "hB" << w << "hL";
216  break;
217  case fileID::positionFile:
218  files(fid) << w << "p" << w << ' ' << w << ' ' << " ";
219  break;
220  }
221  }
222  files(fid).endl();
223 }
224 
225 
226 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
227 
229 (
230  const word& name,
231  const Time& runTime,
232  const dictionary& dict
233 )
234 :
236  logFiles(obr_, name),
237  liquid_(true),
238  alphaName_("alpha"),
239  interpolationScheme_("cellPoint"),
240  direction_(vector::zero),
241  locations_()
242 {
243  read(dict);
244  resetNames({"height", "position"});
245 }
246 
247 
248 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
249 
251 {
252  dict.readIfPresent("alpha", alphaName_);
253  dict.readIfPresent("liquid", liquid_);
254  dict.readEntry("locations", locations_);
255  dict.readIfPresent("interpolationScheme", interpolationScheme_);
256  dict.readIfPresent("direction", direction_);
257 
258  return true;
259 }
260 
261 
263 {
264  return true;
265 }
266 
267 
269 {
270  return true;
271 }
272 
273 
275 {
276  logFiles::write();
277 
278  writePositions();
279 
280  return true;
281 }
282 
283 
284 // ************************************************************************* //
volFields.H
runTime
engineTime & runTime
Definition: createEngineTime.H:13
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::boundBox::mag
scalar mag() const
The magnitude of the bounding box span.
Definition: boundBoxI.H:133
Foam::BitOps::set
void set(List< bool > &bools, const labelRange &range)
Set the specified range 'on' in a boolList.
Definition: BitOps.C:37
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::functionObjects::interfaceHeight::write
virtual bool write()
Write.
Definition: interfaceHeight.C:274
Foam::functionObjects::interfaceHeight::read
virtual bool read(const dictionary &)
Read.
Definition: interfaceHeight.C:250
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::functionObjects::interfaceHeight::execute
virtual bool execute()
Execute.
Definition: interfaceHeight.C:262
Foam::read
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
Foam::functionObjects::logFiles::files
PtrList< OFstream > & files()
Return access to the files.
Definition: logFiles.C:116
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::UPstream::master
static bool master(const label communicator=worldComm)
Am I the master process.
Definition: UPstream.H:457
interpolation.H
Foam::functionObjects::fvMeshFunctionObject
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
Definition: fvMeshFunctionObject.H:64
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Definition: PstreamReduceOps.H:51
Foam::constant::physicoChemical::c1
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
Foam::interpolation< scalar >::New
static autoPtr< interpolation< scalar > > New(const word &interpolationType, const GeometricField< scalar, fvPatchField, volMesh > &psi)
Return a reference to the specified interpolation scheme.
Definition: interpolationNew.C:36
Foam::functionObjects::writeFile::valueWidth
Omanip< int > valueWidth(const label offset=0) const
Return the value width when writing to stream with optional offset.
Definition: writeFile.C:145
Foam::functionObjects::interfaceHeight::end
virtual bool end()
Execute at the final time-loop.
Definition: interfaceHeight.C:268
Foam::uniformDimensionedVectorField
UniformDimensionedField< vector > uniformDimensionedVectorField
Definition: uniformDimensionedFields.H:50
IOmanip.H
Istream and Ostream manipulators taking arguments.
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::objectRegistry::lookupObject
const Type & lookupObject(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:434
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
fvMesh.H
g
const uniformDimensionedVectorField & g
Definition: createFluidFields.H:26
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
interfaceHeight.H
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::functionObjects::logFiles
functionObject base class for creating, maintaining and writing log files e.g. integrated or averaged...
Definition: logFiles.H:59
Foam::polyMesh::bounds
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:450
Time.H
uniformDimensionedFields.H
Foam::constant::atomic::a0
const dimensionedScalar a0
Bohr radius: default SI units: [m].
meshSearch.H
Foam::tab
constexpr char tab
Definition: Ostream.H:403
Foam::functionObjects::addToRunTimeSelectionTable
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
Foam::Omanip
An Ostream manipulator taking arguments.
Definition: IOmanip.H:50
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::functionObjects::fvMeshFunctionObject::mesh_
const fvMesh & mesh_
Reference to the fvMesh.
Definition: fvMeshFunctionObject.H:73
midPointAndFaceSet.H
Foam::functionObjects::defineTypeNameAndDebug
defineTypeNameAndDebug(ObukhovLength, 0)
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
Foam::VectorSpace< Vector< scalar >, scalar, 3 >::zero
static const Vector< scalar > zero
Definition: VectorSpace.H:115
p0
const volScalarField & p0
Definition: EEqn.H:36
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::functionObjects::interfaceHeight::interfaceHeight
interfaceHeight(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: interfaceHeight.C:229
Foam::functionObjects::logFiles::write
virtual bool write()
Write function.
Definition: logFiles.C:149
Foam::functionObjects::interfaceHeight::writeFileHeader
virtual void writeFileHeader(const fileID fid)
Output file header information.
Definition: interfaceHeight.C:155