laminar.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 | Copyright (C) 2011-2016 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 \*---------------------------------------------------------------------------*/
25 
26 #include "laminar.H"
28 #include "fvMesh.H"
29 #include "fvMatrices.H"
30 #include "Time.H"
31 #include "volFields.H"
32 #include "fvmSup.H"
33 #include "kinematicSingleLayer.H"
35 
36 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40 namespace regionModels
41 {
42 namespace surfaceFilmModels
43 {
44 
45 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
46 
47 defineTypeNameAndDebug(laminar, 0);
48 addToRunTimeSelectionTable(filmTurbulenceModel, laminar, dictionary);
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
52 laminar::laminar
53 (
54  surfaceFilmModel& owner,
55  const dictionary& dict
56 )
57 :
58  filmTurbulenceModel(type(), owner, dict),
59  Cf_(readScalar(coeffDict_.lookup("Cf")))
60 {}
61 
62 
63 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
64 
66 {}
67 
68 
69 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
70 
72 {
74  (
75  new volVectorField
76  (
77  IOobject
78  (
79  typeName + ":Us",
84  ),
87  extrapolatedCalculatedFvPatchVectorField::typeName
88  )
89  );
90 
91  // apply quadratic profile
92  tUs.ref() = Foam::sqrt(2.0)*owner_.U();
93  tUs.ref().correctBoundaryConditions();
94 
95  return tUs;
96 }
97 
98 
100 {
101  return tmp<volScalarField>
102  (
103  new volScalarField
104  (
105  IOobject
106  (
107  typeName + ":mut",
109  owner_.regionMesh(),
112  ),
113  owner_.regionMesh(),
115  )
116  );
117 }
118 
119 
121 {
122  // do nothing
123 }
124 
125 
127 {
128  // local reference to film model
129  const kinematicSingleLayer& film =
130  static_cast<const kinematicSingleLayer&>(owner_);
131 
132  // local references to film fields
133  const volScalarField& mu = film.mu();
134  const volVectorField& Uw = film.Uw();
135  const volScalarField& delta = film.delta();
136  const volVectorField& Up = film.UPrimary();
137  const volScalarField& rhop = film.rhoPrimary();
138 
139  // employ simple coeff-based model
140  volScalarField Cs("Cs", Cf_*rhop*mag(Up - U));
141  volScalarField Cw("Cw", mu/(0.3333*(delta + film.deltaSmall())));
142  Cw.min(5000.0);
143 
144  return
145  (
146  - fvm::Sp(Cs, U) + Cs*Up // surface contribution
147  - fvm::Sp(Cw, U) + Cw*Uw // wall contribution
148  );
149 }
150 
151 
152 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
153 
154 } // End namespace surfaceFilmModels
155 } // End namespace regionModels
156 } // End namespace Foam
157 
158 // ************************************************************************* //
scalar delta
dictionary dict
const volVectorField & UPrimary() const
Velocity / [m/s].
U
Definition: pEqn.H:83
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
surfaceFilmModel & owner_
Reference to the owner surface film model.
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual const volVectorField & Uw() const
Return the film wall velocity [m/s].
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
virtual tmp< fvVectorMatrix > Su(volVectorField &U) const
Return the source for the film momentum equation.
Definition: laminar.C:126
Macros for easy insertion into run-time selection tables.
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< volVectorField > Us() const
Return the film surface velocity.
Definition: laminar.C:71
addToRunTimeSelectionTable(surfaceFilmModel, kinematicSingleLayer, mesh)
const volScalarField & delta() const
Return const access to the film thickness / [m].
void min(const dimensioned< Type > &)
static const zero Zero
Definition: zero.H:91
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
const volScalarField & mu() const
Return const access to the dynamic viscosity / [Pa.s].
const volScalarField & rhoPrimary() const
Density / [kg/m3].
const dimensionedScalar mu
Atomic mass unit.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< volScalarField > mut() const
Return the film turbulence viscosity.
Definition: laminar.C:99
virtual void correct()
Correct/update the model.
Definition: laminar.C:120
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
A special matrix type and solver, designed for finite volume solutions of scalar equations.
dimensioned< scalar > mag(const dimensioned< Type > &)
const dimensionedScalar & deltaSmall() const
Return small delta.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:54
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
defineTypeNameAndDebug(kinematicSingleLayer, 0)
virtual const volVectorField & U() const =0
Return the film velocity [m/s].
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
const dimensionSet dimVelocity