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-2018 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 
44 #include "injectionModelList.H"
45 #include "transferModelList.H"
46 #include "forceList.H"
47 #include "filmTurbulenceModel.H"
48 
49 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50 
51 namespace Foam
52 {
53 namespace regionModels
54 {
55 namespace surfaceFilmModels
56 {
57 
58 // Forward class declarations
59 class filmThermoModel;
60 
61 /*---------------------------------------------------------------------------*\
62  Class kinematicSingleLayer Declaration
63 \*---------------------------------------------------------------------------*/
64 
66 :
68 {
69 private:
70 
71  // Private member functions
72 
73  //- Disallow default bitwise copy construct
75 
76  //- Disallow default bitwise assignment
77  void operator=(const kinematicSingleLayer&);
78 
79 
80 protected:
81 
82  // Protected data
83 
84  // Solution parameters
85 
86  //- Momentum predictor
88 
89  //- Number of outer correctors
91 
92  //- Number of PISO-like correctors
93  label nCorr_;
94 
95  //- Number of non-orthogonal correctors
97 
98  //- Cumulative continuity error
99  scalar cumulativeContErr_;
100 
101  //- Small delta
103 
104  //- Film thickness above which Courant number calculation in valid
105  scalar deltaCoLimit_;
106 
107 
108  // Thermo properties
109 
110  // Fields
111 
112  //- Density [kg/m3]
114 
115  //- Dynamic viscosity [Pa.s]
117 
118  //- Surface tension [m/s2]
120 
121 
122  // Fields
123 
124  //- Film thickness [m]
126 
127  //- Film coverage indicator, 1 = covered, 0 = uncovered []
129 
130  //- Velocity - mean [m/s]
132 
133  //- Velocity - surface [m/s]
135 
136  //- Velocity - wall [m/s]
138 
139  //- Film thickness*density (helper field) [kg/m2]
141 
142  //- Mass flux (includes film thickness) [kg.m/s]
144 
145 
146  // Transfer fields
147 
148  //- Film mass available for transfer to the primary region
150 
151  //- Film mass available for transfer to cloud
153 
154  //- Parcel diameters originating from film to cloud
156 
157 
158  // Source term fields
159 
160  // Film region - registered to the film region mesh
161  // Note: need boundary value mapped from primary region, and then
162  // pushed into the patch internal field
163 
164  //- Momementum [kg/m/s2]
166 
167  //- Pressure [Pa]
169 
170  //- Mass [kg/m2/s]
172 
173 
174  // Primary region - registered to the primary region mesh
175  // Internal use only - not read-in
176 
177  //- Momementum [kg/m/s2]
179 
180  //- Pressure [Pa]
182 
183  //- Mass [kg/m2/s]
185 
186 
187  // Fields mapped from primary region - registered to the film region
188  // Note: need both boundary AND patch internal fields to be mapped
189 
190  //- Velocity [m/s]
192 
193  //- Pressure [Pa]
195 
196  //- Density [kg/m3]
198 
199  //- Viscosity [Pa.s]
201 
202 
203  // Sub-models
204 
205  //- Film thermo model
207 
208  //- Available mass for transfer via sub-models
210 
211  //- Cloud injection
213 
214  //- Transfer with the continuous phase
216 
217  //- Turbulence model
219 
220  //- List of film forces
222 
223 
224  // Checks
225 
226  //- Cumulative mass added via sources [kg]
227  scalar addedMassTotal_;
228 
229 
230  // Protected member functions
231 
232  //- Read control parameters from dictionary
233  virtual bool read();
234 
235  //- Correct the thermo fields
236  virtual void correctThermoFields();
237 
238  //- Reset source term fields
239  virtual void resetPrimaryRegionSourceTerms();
240 
241  //- Transfer thermo fields from the primary region to the film region
242  virtual void transferPrimaryRegionThermoFields();
243 
244  //- Transfer source fields from the primary region to the film region
245  virtual void transferPrimaryRegionSourceFields();
246 
247  //- Explicit pressure source contribution
248  virtual tmp<volScalarField> pu();
249 
250  //- Implicit pressure source coefficient
251  virtual tmp<volScalarField> pp();
252 
253  //- Correct film coverage field
254  virtual void correctAlpha();
255 
256  //- Update the film sub-models
257  virtual void updateSubmodels();
258 
259  //- Continuity check
260  virtual void continuityCheck();
261 
262  //- Update film surface velocities
263  virtual void updateSurfaceVelocities();
264 
265  //- Constrain a film region master/slave boundaries of a field to a
266  // given value
267  template<class Type>
268  void constrainFilmField
269  (
270  Type& field,
271  const typename Type::cmptType& value
272  );
273 
274 
275  // Equations
276 
277  //- Solve continuity equation
278  virtual void solveContinuity();
279 
280  //- Solve for film velocity
282  (
283  const volScalarField& pu,
284  const volScalarField& pp
285  );
286 
287  //- Solve coupled velocity-thickness equations
288  virtual void solveThickness
289  (
290  const volScalarField& pu,
291  const volScalarField& pp,
292  const fvVectorMatrix& UEqn
293  );
294 
295 
296 public:
297 
298  //- Runtime type information
299  TypeName("kinematicSingleLayer");
300 
301 
302  // Constructors
303 
304  //- Construct from components
306  (
307  const word& modelType,
308  const fvMesh& mesh,
309  const dimensionedVector& g,
310  const word& regionType,
311  const bool readFields = true
312  );
313 
314 
315  //- Destructor
316  virtual ~kinematicSingleLayer();
317 
318 
319  // Member Functions
320 
321  // Solution parameters
322 
323  //- Courant number evaluation
324  virtual scalar CourantNumber() const;
325 
326  //- Return the momentum predictor
327  inline const Switch& momentumPredictor() const;
328 
329  //- Return the number of outer correctors
330  inline label nOuterCorr() const;
331 
332  //- Return the number of PISO correctors
333  inline label nCorr() const;
334 
335  //- Return the number of non-orthogonal correctors
336  inline label nNonOrthCorr() const;
337 
338  //- Return small delta
339  inline const dimensionedScalar& deltaSmall() const;
340 
341 
342  // Thermo properties
343 
344  //- Return const access to the dynamic viscosity [Pa.s]
345  inline const volScalarField& mu() const;
346 
347  //- Return const access to the surface tension [kg/s2]
348  inline const volScalarField& sigma() const;
349 
350 
351  // Fields
352 
353  //- Return const access to the film thickness [m]
354  inline const volScalarField& delta() const;
355 
356  //- Return the film coverage, 1 = covered, 0 = uncovered []
357  inline const volScalarField& alpha() const;
358 
359  //- Return the film velocity [m/s]
360  virtual const volVectorField& U() const;
361 
362  //- Return the film surface velocity [m/s]
363  virtual const volVectorField& Us() const;
364 
365  //- Return the film wall velocity [m/s]
366  virtual const volVectorField& Uw() const;
367 
368  //- Return the film thickness*density (helper field) [kg/m3]
369  virtual const volScalarField& deltaRho() const;
370 
371  //- Return the film flux [kg.m/s]
372  virtual const surfaceScalarField& phi() const;
373 
374  //- Return the film density [kg/m3]
375  virtual const volScalarField& rho() const;
376 
377  //- Return the film mean temperature [K]
378  virtual const volScalarField& T() const;
379 
380  //- Return the film surface temperature [K]
381  virtual const volScalarField& Ts() const;
382 
383  //- Return the film wall temperature [K]
384  virtual const volScalarField& Tw() const;
385 
386  //- Return the film surface enthalpy [J/kg]
387  virtual const volScalarField& hs() const;
388 
389  //- Return the film specific heat capacity [J/kg/K]
390  virtual const volScalarField& Cp() const;
391 
392  //- Return the film thermal conductivity [W/m/K]
393  virtual const volScalarField& kappa() const;
394 
395 
396  // Transfer fields - to the primary region
397 
398  //- Return mass transfer source - Eulerian phase only
399  virtual tmp<volScalarField> primaryMassTrans() const;
400 
401  //- Return the film mass available for transfer to cloud
402  virtual const volScalarField& cloudMassTrans() const;
403 
404  //- Return the parcel diameters originating from film to cloud
405  virtual const volScalarField& cloudDiameterTrans() const;
406 
407 
408  // External helper functions
409 
410  //- External hook to add sources to the film
411  virtual void addSources
412  (
413  const label patchi, // patchi on primary region
414  const label facei, // facei of patchi
415  const scalar massSource, // [kg]
416  const vector& momentumSource, // [kg.m/s] (tang'l momentum)
417  const scalar pressureSource, // [kg.m/s] (normal momentum)
418  const scalar energySource = 0 // [J]
419  );
420 
421 
422  // Source fields (read/write access)
423 
424  // Primary region
425 
426  //- Momementum [kg/m/s2]
427  inline volVectorField& USpPrimary();
428 
429  //- Pressure [Pa]
430  inline volScalarField& pSpPrimary();
431 
432  //- Mass [kg/m2/s]
433  inline volScalarField& rhoSpPrimary();
434 
435 
436  // Film region
437 
438  //- Momentum [kg/m/s2]
439  inline volVectorField& USp();
440 
441  //- Pressure [Pa]
442  inline volScalarField& pSp();
443 
444  //- Mass [kg/m2/s]
445  inline volScalarField& rhoSp();
446 
447  //- Momentum [kg/m/s2]
448  inline const volVectorField& USp() const;
449 
450  //- Pressure [Pa]
451  inline const volScalarField& pSp() const;
452 
453  //- Mass [kg/m2/s]
454  inline const volScalarField& rhoSp() const;
455 
456 
457  // Fields mapped from primary region
458 
459  //- Velocity [m/s]
460  inline const volVectorField& UPrimary() const;
461 
462  //- Pressure [Pa]
463  inline const volScalarField& pPrimary() const;
464 
465  //- Density [kg/m3]
466  inline const volScalarField& rhoPrimary() const;
467 
468  //- Viscosity [Pa.s]
469  inline const volScalarField& muPrimary() const;
470 
471 
472  // Sub-models
473 
474  //- Film thermo
475  inline const filmThermoModel& filmThermo() const;
476 
477  //- Injection
478  inline injectionModelList& injection();
479 
480  //- Transfer
481  inline transferModelList& transfer();
482 
483  //- Turbulence
484  inline const filmTurbulenceModel& turbulence() const;
485 
486 
487  // Helper functions
488 
489  //- Return the current film mass
490  inline tmp<volScalarField> mass() const;
491 
492  //- Return the change in film mass due to sources/sinks
493  inline tmp<volScalarField> deltaMass() const;
494 
495  //- Return the gravity normal-to-patch component contribution
496  inline tmp<volScalarField> gNorm() const;
497 
498  //- Return the gravity normal-to-patch component contribution
499  // Clipped so that only non-zero if g & nHat_ < 0
500  inline tmp<volScalarField> gNormClipped() const;
501 
502  //- Return the gravity tangential component contributions
503  inline tmp<volVectorField> gTan() const;
504 
505 
506  // Evolution
507 
508  //- Pre-evolve film hook
509  virtual void preEvolveRegion();
510 
511  //- Evolve the film equations
512  virtual void evolveRegion();
513 
514 
515  // Source fields
516 
517  // Mapped into primary region
518 
519  //- Return total mass source - Eulerian phase only
520  virtual tmp<volScalarField::Internal> Srho() const;
521 
522  //- Return mass source for specie i - Eulerian phase only
524  (
525  const label i
526  ) const;
527 
528  //- Return enthalpy source - Eulerian phase only
529  virtual tmp<volScalarField::Internal> Sh() const;
530 
531 
532  // I-O
533 
534  //- Provide some feedback
535  virtual void info();
536 };
537 
538 
539 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
540 
541 } // End namespace surfaceFilmModels
542 } // End namespace regionModels
543 } // End namespace Foam
544 
545 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
546 
547 #ifdef NoRepository
549 #endif
550 
551 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
552 
553 #include "kinematicSingleLayerI.H"
554 
555 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
556 
557 #endif
558 
559 // ************************************************************************* //
virtual bool read()
Read control parameters from dictionary.
Foam::surfaceFields.
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
List container for film injection models.
virtual const volVectorField & U() const
Return the film velocity [m/s].
autoPtr< filmTurbulenceModel > turbulence_
Turbulence model.
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
virtual const volScalarField & deltaRho() const
Return the film thickness*density (helper field) [kg/m3].
virtual const volScalarField & Ts() const
Return the film surface temperature [K].
volScalarField primaryMassTrans_
Film mass available for transfer to the primary region.
Kinematic form of single-cell layer surface film model.
label nOuterCorr() const
Return the number of outer correctors.
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
label nCorr() const
Return the number of PISO correctors.
virtual void transferPrimaryRegionThermoFields()
Transfer thermo fields from the primary region to the film region.
List container for film sources.
Definition: forceList.H:53
virtual const volVectorField & Us() const
Return the film surface velocity [m/s].
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
virtual tmp< volScalarField::Internal > Srho() const
Return total mass source - Eulerian phase only.
tmp< volScalarField > gNorm() const
Return the gravity normal-to-patch component contribution.
scalarField availableMass_
Available mass for transfer via sub-models.
volScalarField cloudDiameterTrans_
Parcel diameters originating from film to cloud.
scalar addedMassTotal_
Cumulative mass added via sources [kg].
const volScalarField & rhoPrimary() const
Density [kg/m3].
virtual tmp< volScalarField > primaryMassTrans() const
Return mass transfer source - Eulerian phase only.
const Switch & momentumPredictor() const
Return the momentum predictor.
virtual void solveThickness(const volScalarField &pu, const volScalarField &pp, const fvVectorMatrix &UEqn)
Solve coupled velocity-thickness equations.
virtual void evolveRegion()
Evolve the film equations.
virtual void solveContinuity()
Solve continuity equation.
tmp< volScalarField > deltaMass() const
Return the change in film mass due to sources/sinks.
virtual tmp< volScalarField > pp()
Implicit pressure source coefficient.
volScalarField deltaRho_
Film thickness*density (helper field) [kg/m2].
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 tmp< fvVectorMatrix > solveMomentum(const volScalarField &pu, const volScalarField &pp)
Solve for film velocity.
virtual const volScalarField & cloudMassTrans() const
Return the film mass available for transfer to cloud.
tmp< volVectorField > gTan() const
Return the gravity tangential component contributions.
tmp< volScalarField > mass() const
Return the current film mass.
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.
tmp< volScalarField > gNormClipped() const
Return the gravity normal-to-patch component contribution.
const dimensionedScalar & deltaSmall() const
Return small delta.
virtual 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
label nNonOrthCorr_
Number of non-orthogonal correctors.
virtual const volScalarField & kappa() const
Return the film thermal conductivity [W/m/K].
virtual tmp< volScalarField > pu()
Explicit pressure source contribution.
const volScalarField & pPrimary() const
Pressure [Pa].
List container for film transfer models.
const dimensionedVector & g() const
Return the accleration due to gravity.
virtual void correctAlpha()
Correct film coverage field.
virtual const volScalarField & hs() const
Return the film surface enthalpy [J/kg].
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.
virtual void updateSubmodels()
Update the film sub-models.
label patchi
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
surfaceScalarField phi_
Mass flux (includes film thickness) [kg.m/s].
const volScalarField & sigma() const
Return const access to the surface tension [kg/s2].
label nNonOrthCorr() const
Return the number of non-orthogonal correctors.
fvVectorMatrix & UEqn
Definition: UEqn.H:13
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered [].
virtual void updateSurfaceVelocities()
Update film surface velocities.
virtual const volScalarField & Tw() const
Return the film wall temperature [K].
A special matrix type and solver, designed for finite volume solutions of scalar equations.
TypeName("kinematicSingleLayer")
Runtime type information.
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.
virtual const volScalarField & rho() const
Return the film density [kg/m3].
volScalarField cloudMassTrans_
Film mass available for transfer to cloud.
void constrainFilmField(Type &field, const typename Type::cmptType &value)
Constrain a film region master/slave boundaries of a field to a.
A class for managing temporary objects.
Definition: PtrList.H:53
virtual const volScalarField & Cp() const
Return the film specific heat capacity [J/kg/K].
const filmTurbulenceModel & turbulence() const
Turbulence.
virtual const volScalarField & T() const
Return the film mean temperature [K].
virtual tmp< volScalarField::Internal > Sh() const
Return enthalpy source - Eulerian phase only.
virtual scalar CourantNumber() const
Courant number evaluation.
scalar deltaCoLimit_
Film thickness above which Courant number calculation in valid.
const volScalarField & alpha() const
Return the film coverage, 1 = covered, 0 = uncovered [].
Namespace for OpenFOAM.
virtual void correctThermoFields()
Correct the thermo fields.