laminar.C
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-2021 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);
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
53 (
55  const dictionary& dict
56 )
57 :
58  momentumTransportModel(type(), film, dict),
59  Cf_(coeffDict_.lookup<scalar>("Cf"))
60 {}
61 
62 
63 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
64 
66 {}
67 
68 
69 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
70 
72 {
73  // Evaluate surface velocity assuming parabolic profile
75  (
77  (
78  IOobject::modelName("Us", typeName),
79  1.5*filmModel_.U()
80  )
81  );
82 
83  return tUs;
84 }
85 
86 
88 {}
89 
90 
92 {
93  // local reference to film model
94  const kinematicSingleLayer& film =
95  static_cast<const kinematicSingleLayer&>(filmModel_);
96 
97  // local references to film fields
98  const volScalarField::Internal& mu = film.mu();
99  const volScalarField::Internal& rho = film.rho();
100  const volScalarField::Internal& delta = film.delta();
101  const volVectorField::Internal& Up = film.UPrimary();
102  const volScalarField::Internal& rhop = film.rhoPrimary();
103  const volScalarField::Internal& VbyA = film.VbyA();
104 
105  // Employ simple coeff-based model
106  volScalarField::Internal Cs("Cs", Cf_*rhop*mag(Up - U)/VbyA);
108  (
109  "Cw",
110  mu/((1.0/3.0)*VbyA*delta + 1e-5*mu*film.time().deltaT()/rho)
111  );
112 
113  return
114  (
115  - fvm::Sp(Cs, U) + Cs*Up // surface contribution
116  - fvm::Sp(Cw, U) + Cw*film.Uw() // wall contribution
117  );
118 }
119 
120 
121 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
122 
123 } // End namespace surfaceFilmModels
124 } // End namespace regionModels
125 } // End namespace Foam
126 
127 // ************************************************************************* //
scalar delta
dictionary dict
const volScalarField & VbyA() const
Return the cell layer volume/area [m].
virtual tmp< volVectorField::Internal > Us() const
Return the film surface velocity.
Definition: laminar.C:71
static tmp< DimensionedField< Type, GeoMesh > > New(const word &name, const Mesh &mesh, const dimensionSet &)
Return a temporary field constructed from name, mesh.
Kinematic form of single-cell layer surface film model.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:37
const volScalarField & rhoPrimary() const
Density [kg/m^3].
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
laminar(surfaceFilmRegionModel &film, const dictionary &dict)
Construct from surface film model.
Definition: laminar.C:53
Macros for easy insertion into run-time selection tables.
CompressibleMomentumTransportModel< dynamicTransportModel > momentumTransportModel
tmp< volVectorField::Internal > Uw() const
Return the film wall velocity [m/s].
const volVectorField & UPrimary() const
Velocity [m/s].
const volScalarField & mu() const
Return const access to the dynamic viscosity [Pa.s].
virtual tmp< fvVectorMatrix > Su(volVectorField &U) const
Return the source for the film momentum equation.
Definition: laminar.C:91
virtual const volVectorField & U() const =0
Return the film velocity [m/s].
const dimensionedScalar mu
Atomic mass unit.
const volScalarField & delta() const
Return const access to the film thickness [m].
U
Definition: pEqn.H:72
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void correct()
Correct/update the model.
Definition: laminar.C:87
A special matrix type and solver, designed for finite volume solutions of scalar equations.
surfaceFilmRegionModel & filmModel_
Reference to the film surface film model.
dimensioned< scalar > mag(const dimensioned< Type > &)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
const volScalarField & rho() const
Return the film density [kg/m^3].
static word modelName(Name name, const word &model)
Return the name of the object within the given model.
A class for managing temporary objects.
Definition: PtrList.H:53
defineTypeNameAndDebug(kinematicSingleLayer, 0)
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
Calculate the matrix for implicit and explicit sources.
addToRunTimeSelectionTable(surfaceFilmModel, noFilm, mesh)
Namespace for OpenFOAM.