kinematicSingleLayerI.H
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) 2011-2017 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 "filmThermoModel.H"
30 #include "surfaceInterpolate.H"
31 #include "fvcSurfaceIntegrate.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace regionModels
38 {
39 namespace surfaceFilmModels
40 {
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
45 {
46  return momentumPredictor_;
47 }
48 
49 
51 {
52  return nOuterCorr_;
53 }
54 
55 
56 inline label kinematicSingleLayer::nCorr() const
57 {
58  return nCorr_;
59 }
60 
61 
63 {
64  return nNonOrthCorr_;
65 }
66 
67 
69 {
70  return deltaSmall_;
71 }
72 
73 
75 {
76  return mu_;
77 }
78 
79 
81 {
82  return sigma_;
83 }
84 
85 
87 {
88  return delta_;
89 }
90 
91 
93 {
94  return alpha_;
95 }
96 
97 
99 {
100  return USpPrimary_;
101 }
102 
103 
105 {
106  return pSpPrimary_;
107 }
108 
109 
111 {
112  return rhoSpPrimary_;
113 }
114 
115 
117 {
118  return USp_;
119 }
120 
121 
123 {
124  return pSp_;
125 }
126 
127 
129 {
130  return rhoSp_;
131 }
132 
133 
135 {
136  return USp_;
137 }
138 
139 
141 {
142  return pSp_;
143 }
144 
145 
147 {
148  return rhoSp_;
149 }
150 
151 
153 {
154  return UPrimary_;
155 }
156 
157 
159 {
160  return pPrimary_;
161 }
162 
163 
165 {
166  return rhoPrimary_;
167 }
168 
169 
171 {
172  return muPrimary_;
173 }
174 
175 
177 {
178  return *filmThermo_;
179 }
180 
181 
183 {
184  return injection_;
185 }
186 
187 
189 {
190  return transfer_;
191 }
192 
193 
195 {
196  return *turbulence_;
197 }
198 
199 
201 {
202  return deltaRho_*magSf();
203 }
204 
205 
207 {
208  return rhoSp_*magSf()*time().deltaT();
209 }
210 
211 
213 {
214  tmp<volScalarField> tgNorm
215  (
216  new volScalarField
217  (
218  IOobject
219  (
220  "gNorm",
221  time().timeName(),
222  regionMesh(),
225  ),
226  g_ & nHat()
227  )
228  );
229 
230  return tgNorm;
231 }
232 
233 
235 {
236  tmp<volScalarField> tgNormClipped
237  (
238  new volScalarField
239  (
240  IOobject
241  (
242  "gNormClipped",
243  time().timeName(),
244  regionMesh(),
247  ),
248  g_ & nHat()
249  )
250  );
251 
252  volScalarField& gNormClipped = tgNormClipped.ref();
253  gNormClipped.min(0.0);
254 
255  return tgNormClipped;
256 }
257 
258 
260 {
261  tmp<volVectorField> tgTan
262  (
263  new volVectorField
264  (
265  IOobject
266  (
267  "gTan",
268  time().timeName(),
269  regionMesh(),
272  ),
273  g_ - nHat()*gNorm()
274  )
275  );
276 
277  return tgTan;
278 }
279 
281 (
282  const label patchI
283 ) const
284 {
285  const vectorField& nH = nHat().boundaryField()[patchI];
286  const vector& g = g_.value();
287  tmp<vectorField> tgTan(new vectorField(g - nH*(g & nH)));
288 
289  return tgTan;
290 }
291 
292 
293 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
294 
295 } // End namespace surfaceFilmModels
296 } // End namespace regionModels
297 } // End namespace Foam
298 
299 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSpPrimary_
volScalarField rhoSpPrimary_
Mass [kg/m2/s].
Definition: kinematicSingleLayer.H:186
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USp
volVectorField & USp()
Momentum [kg/m/s2].
Definition: kinematicSingleLayerI.H:116
Foam::Switch
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:77
filmThermoModel.H
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USpPrimary_
volVectorField USpPrimary_
Momentum [kg/m/s2].
Definition: kinematicSingleLayer.H:180
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USp_
volVectorField USp_
Momentum [kg/m/s2].
Definition: kinematicSingleLayer.H:167
Foam::TimeState::deltaT
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:55
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::turbulence
const filmTurbulenceModel & turbulence() const
Turbulence.
Definition: kinematicSingleLayerI.H:194
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::injection_
injectionModelList injection_
Cloud injection.
Definition: kinematicSingleLayer.H:214
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::turbulence_
autoPtr< filmTurbulenceModel > turbulence_
Turbulence model.
Definition: kinematicSingleLayer.H:220
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::filmThermo_
autoPtr< filmThermoModel > filmThermo_
Film thermo model.
Definition: kinematicSingleLayer.H:208
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaSmall_
const dimensionedScalar deltaSmall_
Small delta.
Definition: kinematicSingleLayer.H:104
Foam::regionModels::surfaceFilmModels::filmTurbulenceModel
Base class for film turbulence models.
Definition: filmTurbulenceModel.H:58
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nCorr
label nCorr() const
Return the number of PISO correctors.
Definition: kinematicSingleLayerI.H:56
Foam::regionModels::regionModel::time
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:40
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::mu
const volScalarField & mu() const
Return const access to the dynamic viscosity [Pa.s].
Definition: kinematicSingleLayerI.H:74
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::sigma_
volScalarField sigma_
Surface tension [m/s2].
Definition: kinematicSingleLayer.H:121
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSp
volScalarField & rhoSp()
Mass [kg/m2/s].
Definition: kinematicSingleLayerI.H:128
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::transfer
transferModelList & transfer()
Transfer.
Definition: kinematicSingleLayerI.H:188
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::UPrimary
const volVectorField & UPrimary() const
Velocity [m/s].
Definition: kinematicSingleLayerI.H:152
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pPrimary_
volScalarField pPrimary_
Pressure [Pa].
Definition: kinematicSingleLayer.H:196
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaSmall
const dimensionedScalar & deltaSmall() const
Return small delta.
Definition: kinematicSingleLayerI.H:68
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::muPrimary_
volScalarField muPrimary_
Viscosity [Pa.s].
Definition: kinematicSingleLayer.H:202
Foam::regionModels::surfaceFilmModels::transferModelList
Definition: transferModelList.H:56
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:54
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::gNormClipped
tmp< volScalarField > gNormClipped() const
Return the gravity normal-to-patch component contribution.
Definition: kinematicSingleLayerI.H:234
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoPrimary
const volScalarField & rhoPrimary() const
Density [kg/m3].
Definition: kinematicSingleLayerI.H:164
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::injection
injectionModelList & injection()
Injection.
Definition: kinematicSingleLayerI.H:182
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaRho_
volScalarField deltaRho_
Film thickness*density (helper field) [kg/m2].
Definition: kinematicSingleLayer.H:142
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::Field< vector >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::gNorm
tmp< volScalarField > gNorm() const
Return the gravity normal-to-patch component contribution.
Definition: kinematicSingleLayerI.H:212
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::delta_
volScalarField delta_
Film thickness [m].
Definition: kinematicSingleLayer.H:127
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::deltaMass
tmp< volScalarField > deltaMass() const
Return the change in film mass due to sources/sinks.
Definition: kinematicSingleLayerI.H:206
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::USpPrimary
volVectorField & USpPrimary()
Momentum [kg/m/s2].
Definition: kinematicSingleLayerI.H:98
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::filmThermo
const filmThermoModel & filmThermo() const
Film thermo.
Definition: kinematicSingleLayerI.H:176
Foam::regionModels::singleLayerRegion::magSf
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
Definition: singleLayerRegion.C:223
Foam::regionModels::regionModel::regionMesh
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:64
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::gTan
tmp< volVectorField > gTan() const
Return the gravity tangential component contributions.
Definition: kinematicSingleLayerI.H:259
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pPrimary
const volScalarField & pPrimary() const
Pressure [Pa].
Definition: kinematicSingleLayerI.H:158
fvcSurfaceIntegrate.H
Surface integrate surfaceField creating a volField. Surface sum a surfaceField creating a volField.
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::muPrimary
const volScalarField & muPrimary() const
Viscosity [Pa.s].
Definition: kinematicSingleLayerI.H:170
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::momentumPredictor_
Switch momentumPredictor_
Momentum predictor.
Definition: kinematicSingleLayer.H:89
Foam::dimensioned< scalar >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nNonOrthCorr
label nNonOrthCorr() const
Return the number of non-orthogonal correctors.
Definition: kinematicSingleLayerI.H:62
g
const uniformDimensionedVectorField & g
Definition: createFluidFields.H:26
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nCorr_
label nCorr_
Number of PISO-like correctors.
Definition: kinematicSingleLayer.H:95
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::mass
tmp< volScalarField > mass() const
Return the current film mass.
Definition: kinematicSingleLayerI.H:200
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::UPrimary_
volVectorField UPrimary_
Velocity [m/s].
Definition: kinematicSingleLayer.H:193
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSp_
volScalarField pSp_
Pressure [Pa].
Definition: kinematicSingleLayer.H:170
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::momentumPredictor
Switch momentumPredictor() const
Return the momentum predictor.
Definition: kinematicSingleLayerI.H:44
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::sigma
const volScalarField & sigma() const
Return const access to the surface tension [kg/s2].
Definition: kinematicSingleLayerI.H:80
Foam::regionModels::surfaceFilmModels::injectionModelList
List container for film injection models.
Definition: injectionModelList.H:56
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSpPrimary
volScalarField & rhoSpPrimary()
Mass [kg/m2/s].
Definition: kinematicSingleLayerI.H:110
Foam::regionModels::surfaceFilmModels::filmThermoModel
Base class for film thermo models.
Definition: filmThermoModel.H:56
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSp
volScalarField & pSp()
Pressure [Pa].
Definition: kinematicSingleLayerI.H:122
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::alpha_
volScalarField alpha_
Film coverage indicator, 1 = covered, 0 = uncovered [].
Definition: kinematicSingleLayer.H:130
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSpPrimary_
volScalarField pSpPrimary_
Pressure [Pa].
Definition: kinematicSingleLayer.H:183
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::mu_
volScalarField mu_
Dynamic viscosity [Pa.s].
Definition: kinematicSingleLayer.H:118
Foam::Vector< scalar >
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nNonOrthCorr_
label nNonOrthCorr_
Number of non-orthogonal correctors.
Definition: kinematicSingleLayer.H:98
Foam::regionModels::singleLayerRegion::nHat
virtual const volVectorField & nHat() const
Return the patch normal vectors.
Definition: singleLayerRegion.C:210
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoPrimary_
volScalarField rhoPrimary_
Density [kg/m3].
Definition: kinematicSingleLayer.H:199
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::delta
const volScalarField & delta() const
Return const access to the film thickness [m].
Definition: kinematicSingleLayerI.H:86
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::alpha
const volScalarField & alpha() const
Return the film coverage, 1 = covered, 0 = uncovered [].
Definition: kinematicSingleLayerI.H:92
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::transfer_
transferModelList transfer_
Transfer with the continuous phase.
Definition: kinematicSingleLayer.H:217
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::regionModels::surfaceFilmModels::surfaceFilmRegionModel::g_
const dimensionedVector & g_
Acceleration due to gravity [m/s2].
Definition: surfaceFilmRegionModel.H:74
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::rhoSp_
volScalarField rhoSp_
Mass [kg/m2/s].
Definition: kinematicSingleLayer.H:173
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::pSpPrimary
volScalarField & pSpPrimary()
Pressure [Pa].
Definition: kinematicSingleLayerI.H:104
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nOuterCorr_
label nOuterCorr_
Number of outer correctors.
Definition: kinematicSingleLayer.H:92
Foam::regionModels::surfaceFilmModels::kinematicSingleLayer::nOuterCorr
label nOuterCorr() const
Return the number of outer correctors.
Definition: kinematicSingleLayerI.H:50