liquidFilmBase.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) 2020-2023 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "liquidFilmBase.H"
29 #include "gravityMeshObject.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace regionModels
39 {
40 namespace areaSurfaceFilmModels
41 {
42 
43 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44 
45 defineTypeNameAndDebug(liquidFilmBase, 0);
46 
48 
49 const Foam::word liquidFilmName("liquidFilm");
50 
51 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52 
54 (
55  const word& modelType,
56  const fvMesh& mesh,
57  const dictionary& dict
58 )
59 :
60  regionFaModel(mesh, liquidFilmName, modelType, dict, true),
61  momentumPredictor_
62  (
63  this->solution().subDict("PIMPLE").get<bool>("momentumPredictor")
64  ),
65  nOuterCorr_
66  (
67  this->solution().subDict("PIMPLE").get<label>("nOuterCorr")
68  ),
69  nCorr_(this->solution().subDict("PIMPLE").get<label>("nCorr")),
70  nFilmCorr_
71  (
72  this->solution().subDict("PIMPLE").get<label>("nFilmCorr")
73  ),
74  h0_("h0", dimLength, 1e-7, dict),
75  deltaWet_("deltaWet", dimLength, 1e-4, dict),
76  UName_(dict.get<word>("U")),
77  pName_(dict.getOrDefault<word>("p", word::null)),
78  pRef_(dict.get<scalar>("pRef")),
79  h_
80  (
81  IOobject
82  (
83  "hf_" + regionName_,
84  regionMesh().time().timeName(),
85  regionMesh().thisDb(),
86  IOobject::MUST_READ,
87  IOobject::AUTO_WRITE
88  ),
89  regionMesh()
90  ),
91  Uf_
92  (
93  IOobject
94  (
95  "Uf_" + regionName_,
96  regionMesh().time().timeName(),
97  regionMesh().thisDb(),
98  IOobject::MUST_READ,
99  IOobject::AUTO_WRITE
100  ),
101  regionMesh()
102  ),
103  pf_
104  (
105  IOobject
106  (
107  "pf_" + regionName_,
108  regionMesh().time().timeName(),
109  regionMesh().thisDb(),
110  IOobject::NO_READ,
111  IOobject::AUTO_WRITE
112  ),
113  regionMesh(),
115  ),
116  ppf_
117  (
118  IOobject
119  (
120  "ppf_" + regionName_,
121  regionMesh().time().timeName(),
122  regionMesh().thisDb(),
123  IOobject::NO_READ,
124  IOobject::NO_WRITE
125  ),
126  regionMesh(),
128  ),
129  phif_
130  (
131  IOobject
132  (
133  "phif_" + regionName_,
134  regionMesh().time().timeName(),
135  regionMesh().thisDb(),
136  IOobject::READ_IF_PRESENT,
137  IOobject::AUTO_WRITE
138  ),
139  fac::interpolate(Uf_) & regionMesh().Le()
140  ),
141  phi2s_
142  (
143  IOobject
144  (
145  "phi2s_" + regionName_,
146  regionMesh().time().timeName(),
147  regionMesh().thisDb(),
148  IOobject::READ_IF_PRESENT,
149  IOobject::AUTO_WRITE
150  ),
151  fac::interpolate(h_*Uf_) & regionMesh().Le()
152  ),
153  gn_
154  (
155  IOobject
156  (
157  "gn",
158  regionMesh().time().timeName(),
159  regionMesh().thisDb(),
160  IOobject::NO_READ,
161  IOobject::NO_WRITE
162  ),
163  regionMesh(),
165  ),
166  g_(meshObjects::gravity::New(primaryMesh().time())),
167  massSource_
168  (
169  IOobject
170  (
171  "massSource",
172  primaryMesh().time().timeName(),
173  primaryMesh().thisDb()
174  ),
175  primaryMesh(),
177  ),
178  momentumSource_
179  (
180  IOobject
181  (
182  "momentumSource",
183  primaryMesh().time().timeName(),
184  primaryMesh().thisDb()
185  ),
186  primaryMesh(),
188  ),
189  pnSource_
190  (
191  IOobject
192  (
193  "pnSource",
194  primaryMesh().time().timeName(),
195  primaryMesh().thisDb()
196  ),
197  primaryMesh(),
199  ),
200  addedMassTotal_(0),
201  faOptions_(Foam::fa::options::New(primaryMesh()))
202 {
204 
205  gn_ = g_ & ns;
206 
207  if (!faOptions_.optionList::size())
208  {
209  Info << "No finite area options present" << endl;
210  }
211 }
212 
213 
214 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
217 {}
218 
219 
220 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
221 
222 scalar liquidFilmBase::CourantNumber() const
223 {
224  scalar CoNum = 0.0;
225  scalar velMag = 0.0;
226 
227  edgeScalarField SfUfbyDelta
228  (
230  );
231 
232  CoNum =
233  max(SfUfbyDelta/regionMesh().magLe()).value()*time().deltaTValue();
234 
235  velMag = max(mag(phif_)/regionMesh().magLe()).value();
236 
239 
240  Info<< "Max film Courant Number: " << CoNum
241  << " Film velocity magnitude: " << velMag << endl;
242 
243  return CoNum;
244 }
245 
246 
248 {
249  auto tUw = areaVectorField::New
250  (
251  "tUw",
253  regionMesh(),
255  );
256  auto& Uw = tUw.ref();
257 
258  if (primaryMesh().moving())
259  {
261 
262  // Set up mapping values per patch
263 
264  PtrMap<vectorField> patchValues(2*patches.size());
265 
266  for (const label patchi : patches)
267  {
268  const auto* wpp = isA<movingWallVelocityFvPatchVectorField>
269  (
270  primaryMesh().boundaryMesh()[patchi]
271  );
272 
273  if (wpp)
274  {
275  patchValues.set(patchi, wpp->Uwall());
276  }
277  }
278 
279  if (patchValues.size())
280  {
281  // Map Up to area
282  tmp<vectorField> UsWall = vsmPtr_->mapToSurface(patchValues);
283 
284  const vectorField& nHat =
286 
287  Uw.primitiveFieldRef() = UsWall() - nHat*(UsWall() & nHat);
288  }
289  }
290 
291  return tUw;
292 }
293 
294 
296 {
297  auto tUs = areaVectorField::New
298  (
299  "tUs",
301  regionMesh(),
303  );
304 
305  // Uf quadratic profile
306  tUs.ref() = Foam::sqrt(2.0)*Uf_;
307 
308  return tUs;
309 }
310 
311 
313 {
314  const volVectorField& Uprimary =
316 
317  auto tUp = areaVectorField::New
318  (
319  "tUp",
321  regionMesh(),
323  );
324  auto& Up = tUp.ref();
325 
326 
327  // Set up mapping values per patch
328 
330 
331  PtrMap<vectorField> patchValues(2*patches.size());
332 
333  // U tangential on primary
334  for (const label patchi : patches)
335  {
336  const fvPatchVectorField& Uw = Uprimary.boundaryField()[patchi];
337 
338  patchValues.set(patchi, -Uw.snGrad());
339  }
340 
341 
342  // Map U tang to surface
343  vsmPtr_->mapToSurface(patchValues, Up.primitiveFieldRef());
344 
345  Up.primitiveFieldRef() *= h_.primitiveField();
346 
347  // U tangent on surface
350  Up.primitiveFieldRef() -= nHat*(Up.primitiveField() & nHat);
351 
352  return tUp;
353 }
354 
355 
357 {
358  auto tpg = areaScalarField::New
359  (
360  "tpg",
362  regionMesh(),
364  );
365  auto& pfg = tpg.ref();
366 
367  if (!pName_.empty())
368  {
369  vsmPtr_->mapInternalToSurface
370  (
371  primaryMesh().lookupObject<volScalarField>(pName_),
372  pfg.primitiveFieldRef()
373  );
374  }
375 
376  return tpg;
377 }
378 
379 
381 {
382  auto talpha = areaScalarField::New
383  (
384  "talpha",
386  regionMesh(),
388  );
389  auto& alpha = talpha.ref();
390 
392 
393  return talpha;
394 }
395 
396 
398 (
399  const label patchi,
400  const label facei,
401  const scalar massSource,
402  const vector& momentumSource,
403  const scalar pressureSource,
404  const scalar energySource
405 )
406 {
407  massSource_.boundaryFieldRef()[patchi][facei] += massSource;
408  pnSource_.boundaryFieldRef()[patchi][facei] += pressureSource;
409  momentumSource_.boundaryFieldRef()[patchi][facei] += momentumSource;
410 }
411 
414 {
416 }
417 
418 
420 {
421  if (debug && primaryMesh().time().writeTime())
422  {
423  massSource_.write();
424  pnSource_.write();
426  }
427 
431 
433 }
434 
435 
436 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
437 
438 } // End namespace areaSurfaceFilmModels
439 } // End namespace regionModels
440 } // End namespace Foam
441 
442 // ************************************************************************* //
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
dictionary dict
const edgeScalarField & deltaCoeffs() const
Return reference to difference factors array.
const Type & lookupObject(const word &name, const bool recursive=false) const
Lookup and return const reference to the object of the given Type. Fatal if not found or the wrong ty...
Calculate the second temporal derivative.
const Internal::FieldType & primitiveField() const noexcept
Return a const-reference to the internal field values.
scalar deltaTValue() const noexcept
Return time step value.
Definition: TimeStateI.H:49
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
tmp< areaVectorField > Up() const
Primary region velocity at film hight. Assume the film to be.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:40
const Internal & internalField() const noexcept
Return a const-reference to the dimensioned internal field.
virtual void postEvolveRegion()
Post-evolve region.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:531
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
A HashTable of pointers to objects of type <T> with a label key.
Definition: HashTableFwd.H:42
const labelList & whichPolyPatches() const
The polyPatches related to the areaMesh, in sorted order.
Definition: faMeshI.H:128
tmp< areaScalarField > pg() const
Map primary static pressure.
const dimensionSet dimless
Dimensionless.
virtual scalar CourantNumber() const
Courant number evaluation.
defineTypeNameAndDebug(kinematicThinFilm, 0)
liquidFilmBase(const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from type name and mesh and dict.
volScalarField pnSource_
Normal pressure by particles.
word timeName
Definition: getTimeIndex.H:3
const dimensionSet dimAcceleration
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
dynamicFvMesh & mesh
defineRunTimeSelectionTable(liquidFilmBase, dictionary)
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
Definition: polyMesh.H:609
A class for handling words, derived from Foam::string.
Definition: word.H:63
label size() const noexcept
The number of entries in the list.
Definition: UPtrListI.H:106
bool set(const Key &key, T *ptr)
Assign a new entry, overwrites existing.
const fvMesh & primaryMesh() const noexcept
Return the reference to the primary mesh database.
virtual bool write(const bool writeOnProc=true) const
Write using setting from DB.
const faMesh & regionMesh() const
Return the region mesh database.
const dimensionSet dimPressure
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
Definition: PtrList.H:159
tmp< areaVectorField > Us() const
Film surface film velocity.
static tmp< GeometricField< vector, faPatchField, areaMesh > > New(const word &name, IOobjectOption::registerOption regOpt, const Mesh &mesh, const dimensionSet &dims, const word &patchFieldType=faPatchField< vector >::calculatedType())
Return tmp field (NO_READ, NO_WRITE) from name, mesh, dimensions and patch type. [Takes current timeN...
scalar CoNum
virtual void preEvolveRegion()
Pre-evolve region.
dimensionedScalar pos0(const dimensionedScalar &ds)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
Base class for area region models.
const Foam::word liquidFilmName("liquidFilm")
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const polyBoundaryMesh & patches
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensionedScalar deltaWet_
Film thickness beyond which face is assumed to be wet.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
const Time & time() const noexcept
Return the reference to the time database.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:92
tmp< areaVectorField > Uw() const
Wall velocity.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: HashPtrTable.H:50
tmp< areaScalarField > alpha() const
Wet indicator using h0.
scalar velMag
Defines the attributes of an object for which implicit objectRegistry management is supported...
Definition: IOobject.H:180
const areaVectorField & faceAreaNormals() const
Return face area normals.
Definition: faMesh.C:1183
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
Do not request registration (bool: false)
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.
autoPtr< volSurfaceMapping > vsmPtr_
Volume/surface mapping.
Namespace for OpenFOAM.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:127
const dimensionSet dimVelocity
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource=0)
Add sources.