All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
kinematicSingleLayer.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::regionModels::surfaceFilmModels::kinematicSingleLayer
26 
27 Description
28  Kinematic form of single-cell layer surface film model
29 
30 SourceFiles
31  kinematicSingleLayer.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef kinematicSingleLayer_H
36 #define kinematicSingleLayer_H
37 
38 #include "surfaceFilmRegionModel.H"
39 #include "fvMesh.H"
40 #include "volFields.H"
41 #include "surfaceFields.H"
42 #include "fvMatrices.H"
43 #include "pimpleControl.H"
44 
45 #include "injectionModelList.H"
46 #include "transferModelList.H"
47 #include "forceList.H"
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 namespace regionModels
55 {
56 namespace surfaceFilmModels
57 {
58 
59 // Forward class declarations
60 class filmThermoModel;
61 
62 /*---------------------------------------------------------------------------*\
63  Class kinematicSingleLayer Declaration
64 \*---------------------------------------------------------------------------*/
65 
67 :
69 {
70 
71 protected:
72 
73  // Protected data
74 
75  // Solution parameters
76 
78 
79  //- Cumulative continuity error
80  scalar cumulativeContErr_;
81 
82  //- Small delta
84 
85  //- Film thickness above which Courant number calculation in valid
86  scalar deltaCoLimit_;
87 
88 
89  // Thermo properties
90 
91  // Fields
92 
93  //- Density [kg/m^3]
95 
96  //- Dynamic viscosity [Pa.s]
98 
99  //- Surface tension [m/s^2]
101 
102 
103  // Fields
104 
105  //- Film thickness [m]
107 
108  //- Film volume fraction in the cell layer []
110 
111  //- Velocity - mean [m/s]
113 
114  //- Velocity - surface [m/s]
116 
117  //- Velocity - wall [m/s]
119 
120  //- Mass flux [kg m/s]
122 
123  //- Film velocity flux [m^3/s]
125 
126  //- Current continuity error caused by delta_ bounding
128 
129  //- Film coverage indicator, 1 = covered, 0 = uncovered []
131 
132 
133  // Transfer fields
134 
135  //- Film mass available for transfer to the primary region
137 
138  //- Film mass available for transfer to cloud
140 
141  //- Parcel diameters originating from film to cloud
143 
144 
145  // Source term fields
146 
147  // Film region - registered to the film region mesh
148  // Note: need boundary value mapped from primary region, and then
149  // pushed into the patch internal field
150 
151  //- Mass [kg/m^2/s]
153 
154  //- Momentum [kg/m/s^2]
156 
157  //- Pressure [Pa]
159 
160 
161  // Primary region - registered to the primary region mesh
162  // Internal use only - not read-in
163 
164  //- Primary region mass source [kg]
166 
167  //- Primary region tangential momentum source [kg m/s]
169 
170  //- Primary region normal momentum source (pressure) [kg m/s]
172 
173 
174  // Fields mapped from primary region - registered to the film region
175  // Note: need both boundary AND patch internal fields to be mapped
176 
177  //- Velocity [m/s]
179 
180  //- Pressure [Pa]
182 
183  //- Density [kg/m^3]
185 
186  //- Viscosity [Pa.s]
188 
189 
190  // Sub-models
191 
192  //- Film thermo model
194 
195  //- Available mass for transfer via sub-models
197 
198  //- Cloud injection
200 
201  //- Transfer with the continuous phase
203 
204  //- Turbulence model
206 
207  //- List of film forces
209 
210 
211  // Checks
212 
213  //- Cumulative mass added via sources [kg]
214  scalar addedMassTotal_;
215 
216 
217  // Protected member functions
218 
219  //- Read control parameters from dictionary
220  virtual bool read();
221 
222  //- Correct the thermo fields
223  virtual void correctThermoFields();
224 
225  //- Reset source term fields
226  virtual void resetPrimaryRegionSourceTerms();
227 
228  //- Transfer thermo fields from the primary region to the film region
229  virtual void transferPrimaryRegionThermoFields();
230 
231  //- Transfer source fields from the primary region to the film region
232  virtual void transferPrimaryRegionSourceFields();
233 
234  //- Hydrostatic pressure coefficient
236 
237  //- Hydrostatic pressure coefficient gradient
239 
240  //- Capillary pressure
242 
243  //- Explicit pressure
245 
246  //- Correct film coverage field
247  virtual void correctCoverage();
248 
249  //- Update the film sub-models
250  virtual void updateSubmodels();
251 
252  // Update continuity error
253  virtual void updateContinuityErr();
254 
255  //- Continuity check
256  virtual void continuityCheck();
257 
258  //- Update film surface velocities
259  virtual void updateSurfaceVelocities();
260 
261  //- Constrain a film region master/slave boundaries of a field to a
262  // given value
263  template<class Type>
265  (
266  const tmp<Type>& field,
267  const typename Type::cmptType& value
268  );
269 
270 
271  // Equations
272 
273  //- Predict delta_ from the continuity equation
274  virtual void predictDelta();
275 
276  //- Solve for film velocity
278  (
279  const volScalarField& pc,
280  const volScalarField& pe
281  );
282 
283  //- Solve for film volume fraction and thickness
284  virtual void solveAlpha
285  (
286  const fvVectorMatrix& UEqn,
287  const volScalarField& pc,
288  const volScalarField& pe
289  );
290 
291 
292 public:
293 
294  //- Runtime type information
295  TypeName("kinematicSingleLayer");
296 
297 
298  // Constructors
299 
300  //- Construct from components
302  (
303  const word& modelType,
304  const fvMesh& mesh,
305  const dimensionedVector& g,
306  const word& regionType,
307  const bool readFields = true
308  );
309 
310  //- Disallow default bitwise copy construction
312 
313 
314  //- Destructor
315  virtual ~kinematicSingleLayer();
316 
317 
318  // Member Functions
319 
320  // Solution parameters
321 
322  //- Courant number evaluation
323  virtual scalar CourantNumber() const;
324 
325  //- Return small delta
326  inline const dimensionedScalar& deltaSmall() const;
327 
328 
329  // Thermo properties
330 
331  //- Return const access to the dynamic viscosity [Pa.s]
332  inline const volScalarField& mu() const;
333 
334  //- Return const access to the surface tension [kg/s^2]
335  inline const volScalarField& sigma() const;
336 
337 
338  // Fields
339 
340  //- Return const access to the film thickness [m]
341  inline const volScalarField& delta() const;
342 
343  //- Return const access to the film volume fraction []
344  inline const volScalarField& alpha() const;
345 
346  //- Return the film density [kg/m^3]
347  inline const volScalarField& rho() const;
348 
349  //- Return the film velocity [m/s]
350  inline const volVectorField& U() const;
351 
352  //- Return the film surface velocity [m/s]
353  inline const volVectorField::Internal& Us() const;
354 
355  //- Return the film wall velocity [m/s]
356  inline const volVectorField::Internal& Uw() const;
357 
358  //- Return the film flux [kg m/s]
359  inline const surfaceScalarField& phi() const;
360 
361  //- Return the film velocity flux [m^3/s]
362  inline const surfaceScalarField& phiU() const;
363 
364  //- Return the current continuity error
365  inline const volScalarField::Internal& continuityErr() const;
366 
367  //- Return the film coverage, 1 = covered, 0 = uncovered []
368  inline const volScalarField& coverage() const;
369 
370 
371  // Transfer fields - to the primary region
372 
373  //- Return mass transfer source - Eulerian phase only
374  virtual tmp<volScalarField> primaryMassTrans() const;
375 
376  //- Return the film mass available for transfer to cloud
377  virtual const volScalarField& cloudMassTrans() const;
378 
379  //- Return the parcel diameters originating from film to cloud
380  virtual const volScalarField& cloudDiameterTrans() const;
381 
382 
383  // External helper functions
384 
385  //- External hook to add sources to the film
386  virtual void addSources
387  (
388  const label patchi, // patchi on primary region
389  const label facei, // facei of patchi
390  const scalar massSource, // [kg]
391  const vector& momentumSource, // [kg m/s] (tang'l momentum)
392  const scalar pressureSource, // [kg m/s] (normal momentum)
393  const scalar energySource = 0 // [J]
394  );
395 
396 
397  // Source fields (read/write access)
398 
399  // Primary region
400 
401  //- Momentum [kg/m/s^2]
402  inline volVectorField& USpPrimary();
403 
404  //- Pressure [Pa]
405  inline volScalarField& pSpPrimary();
406 
407  //- Mass [kg/m^2/s]
408  inline volScalarField& rhoSpPrimary();
409 
410 
411  // Film region
412 
413  //- Momentum [kg/m/s^2]
414  inline volVectorField::Internal& USp();
415 
416  //- Pressure [Pa]
417  inline volScalarField::Internal& pSp();
418 
419  //- Mass [kg/m^2/s]
421 
422  //- Momentum [kg/m/s^2]
423  inline const volVectorField::Internal& USp() const;
424 
425  //- Pressure [Pa]
426  inline const volScalarField::Internal& pSp() const;
427 
428  //- Mass [kg/m^2/s]
429  inline const volScalarField::Internal& rhoSp() const;
430 
431 
432  // Fields mapped from primary region
433 
434  //- Velocity [m/s]
435  inline const volVectorField& UPrimary() const;
436 
437  //- Pressure [Pa]
438  inline const volScalarField& pPrimary() const;
439 
440  //- Density [kg/m^3]
441  inline const volScalarField& rhoPrimary() const;
442 
443  //- Viscosity [Pa.s]
444  inline const volScalarField& muPrimary() const;
445 
446 
447  // Sub-models
448 
449  //- Film thermo
450  inline const filmThermoModel& filmThermo() const;
451 
452  //- Injection
453  inline injectionModelList& injection();
454 
455  //- Transfer
456  inline transferModelList& transfer();
457 
458  //- Turbulence
459  inline const filmMomentumTransportModel& turbulence() const;
460 
461 
462  // Helper functions
463 
464  //- Return the current film mass
465  inline tmp<volScalarField::Internal> mass() const;
466 
467  //- Return the change in film mass due to sources/sinks
469 
470 
471  // Evolution
472 
473  //- Pre-evolve film hook
474  virtual void preEvolveRegion();
475 
476  //- Evolve the film equations
477  virtual void evolveRegion();
478 
479 
480  // Source fields
481 
482  // Mapped into primary region
483 
484  //- Return total mass source - Eulerian phase only
485  virtual tmp<volScalarField::Internal> Srho() const;
486 
487  //- Return mass source for specie i - Eulerian phase only
489  (
490  const label i
491  ) const;
492 
493  //- Return enthalpy source - Eulerian phase only
494  virtual tmp<volScalarField::Internal> Sh() const;
495 
496 
497  // I-O
498 
499  //- Provide some feedback
500  virtual void info();
501 
502 
503  // Member Operators
504 
505  //- Disallow default bitwise assignment
506  void operator=(const kinematicSingleLayer&) = delete;
507 };
508 
509 
510 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
511 
512 } // End namespace surfaceFilmModels
513 } // End namespace regionModels
514 } // End namespace Foam
515 
516 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
517 
518 #ifdef NoRepository
520 #endif
521 
522 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
523 
524 #include "kinematicSingleLayerI.H"
525 
526 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
527 
528 #endif
529 
530 // ************************************************************************* //
virtual bool read()
Read control parameters from dictionary.
Foam::surfaceFields.
const filmMomentumTransportModel & turbulence() const
Turbulence.
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
List container for film injection models.
tmp< volScalarField::Internal > mass() const
Return the current film mass.
const volVectorField::Internal & Uw() const
Return the film wall velocity [m/s].
const volVectorField & U() const
Return the film velocity [m/s].
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
volScalarField primaryMassTrans_
Film mass available for transfer to the primary region.
Kinematic form of single-cell layer surface film model.
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
volScalarField::Internal continuityErr_
Current continuity error caused by delta_ bounding.
void operator=(const kinematicSingleLayer &)=delete
Disallow default bitwise assignment.
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
List container for film sources.
Definition: forceList.H:53
virtual void correctCoverage()
Correct film coverage field.
volVectorField::Internal USp_
Momentum [kg/m/s^2].
volVectorField::Internal Us_
Velocity - surface [m/s].
virtual tmp< volScalarField::Internal > Srho() const
Return total mass source - Eulerian phase only.
scalarField availableMass_
Available mass for transfer via sub-models.
volScalarField cloudDiameterTrans_
Parcel diameters originating from film to cloud.
volScalarField coverage_
Film coverage indicator, 1 = covered, 0 = uncovered [].
scalar addedMassTotal_
Cumulative mass added via sources [kg].
volVectorField USpPrimary_
Primary region tangential momentum source [kg m/s].
const volScalarField & rhoPrimary() const
Density [kg/m^3].
volScalarField::Internal & rhoSp()
Mass [kg/m^2/s].
virtual tmp< volScalarField > primaryMassTrans() const
Return mass transfer source - Eulerian phase only.
tmp< surfaceScalarField > rhog() const
Hydrostatic pressure coefficient.
virtual void evolveRegion()
Evolve the film equations.
virtual tmp< fvVectorMatrix > solveMomentum(const volScalarField &pc, const volScalarField &pe)
Solve for film velocity.
const volScalarField & muPrimary() const
Viscosity [Pa.s].
virtual void transferPrimaryRegionSourceFields()
Transfer source fields from the primary region to the film region.
const filmThermoModel & filmThermo() const
Film thermo.
const volVectorField & UPrimary() const
Velocity [m/s].
dynamicFvMesh & mesh
const volScalarField & mu() const
Return const access to the dynamic viscosity [Pa.s].
virtual const volScalarField & cloudMassTrans() const
Return the film mass available for transfer to cloud.
const surfaceScalarField & phiU() const
Return the film velocity flux [m^3/s].
A class for handling words, derived from string.
Definition: word.H:59
virtual const volScalarField & cloudDiameterTrans() const
Return the parcel diameters originating from film to cloud.
const dimensionedScalar & deltaSmall() const
Return small delta.
const volScalarField & coverage() const
Return the film coverage, 1 = covered, 0 = uncovered [].
const surfaceScalarField & phi() const
Return the film flux [kg m/s].
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
virtual void solveAlpha(const fvVectorMatrix &UEqn, const volScalarField &pc, const volScalarField &pe)
Solve for film volume fraction and thickness.
const volScalarField::Internal & continuityErr() const
Return the current continuity error.
volScalarField pSpPrimary_
Primary region normal momentum source (pressure) [kg m/s].
volVectorField::Internal & USp()
Momentum [kg/m/s^2].
const volScalarField & pPrimary() const
Pressure [Pa].
List container for film transfer models.
const dimensionedVector & g() const
Return the acceleration due to gravity.
virtual void resetPrimaryRegionSourceTerms()
Reset source term fields.
const volScalarField & delta() const
Return const access to the film thickness [m].
transferModelList transfer_
Transfer with the continuous phase.
kinematicSingleLayer(const word &modelType, const fvMesh &mesh, const dimensionedVector &g, const word &regionType, const bool readFields=true)
Construct from components.
virtual void updateSubmodels()
Update the film sub-models.
label patchi
virtual void predictDelta()
Predict delta_ from the continuity equation.
const volVectorField::Internal & Us() const
Return the film surface velocity [m/s].
fvVectorMatrix & UEqn
Definition: UEqn.H:11
Pimple control class. Provides time-loop control methods which exit the simulation once convergence c...
Definition: pimpleControl.H:74
const volScalarField & sigma() const
Return const access to the surface tension [kg/s^2].
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
volScalarField alpha_
Film volume fraction in the cell layer [].
virtual void updateSurfaceVelocities()
Update film surface velocities.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A special matrix type and solver, designed for finite volume solutions of scalar equations.
TypeName("kinematicSingleLayer")
Runtime type information.
tmp< surfaceScalarField > gGradRho() const
Hydrostatic pressure coefficient gradient.
volScalarField rhoSpPrimary_
Primary region mass source [kg].
surfaceScalarField phiU_
Film velocity flux [m^3/s].
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
virtual void addSources(const label patchi, const label facei, const scalar massSource, const vector &momentumSource, const scalar pressureSource, const scalar energySource=0)
External hook to add sources to the film.
const volScalarField & rho() const
Return the film density [kg/m^3].
volVectorField::Internal Uw_
Velocity - wall [m/s].
volScalarField cloudMassTrans_
Film mass available for transfer to cloud.
rDeltaTY field()
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< Type > constrainFilmField(const tmp< Type > &field, const typename Type::cmptType &value)
Constrain a film region master/slave boundaries of a field to a.
virtual tmp< volScalarField::Internal > Sh() const
Return enthalpy source - Eulerian phase only.
virtual scalar CourantNumber() const
Courant number evaluation.
autoPtr< filmMomentumTransportModel > turbulence_
Turbulence model.
scalar deltaCoLimit_
Film thickness above which Courant number calculation in valid.
const volScalarField & alpha() const
Return const access to the film volume fraction [].
Namespace for OpenFOAM.
virtual void correctThermoFields()
Correct the thermo fields.
tmp< volScalarField::Internal > deltaMass() const
Return the change in film mass due to sources/sinks.