waxSolventViscosity.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) 2017-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 \*---------------------------------------------------------------------------*/
25 
26 #include "waxSolventViscosity.H"
27 #include "kinematicSingleLayer.H"
28 #include "waxSolventEvaporation.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace regionModels
36 {
37 namespace surfaceFilmModels
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 defineTypeNameAndDebug(waxSolventViscosity, 0);
43 
45 (
46  filmViscosityModel,
47  waxSolventViscosity,
48  dictionary
49 );
50 
51 
52 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
53 
54 void waxSolventViscosity::correctMu()
55 {
56  const kinematicSingleLayer& film = filmType<kinematicSingleLayer>();
57 
59  (
60  film.regionMesh().lookupObject<uniformDimensionedScalarField>
61  (
62  waxSolventEvaporation::typeName + ":Wwax"
63  )
64  );
65 
66  const uniformDimensionedScalarField Wsolvent
67  (
68  film.regionMesh().lookupObject<uniformDimensionedScalarField>
69  (
70  waxSolventEvaporation::typeName + ":Wsolvent"
71  )
72  );
73 
74  const uniformDimensionedScalarField Ysolvent0
75  (
76  film.regionMesh().lookupObject<uniformDimensionedScalarField>
77  (
78  waxSolventEvaporation::typeName + ":Ysolvent0"
79  )
80  );
81 
82  const volScalarField& Ysolvent
83  (
84  film.regionMesh().lookupObject<volScalarField>
85  (
86  waxSolventEvaporation::typeName + ":Ysolvent"
87  )
88  );
89 
90  const volScalarField Xsolvent
91  (
92  Ysolvent*Wsolvent/((1 - Ysolvent)*Wwax + Ysolvent*Wsolvent)
93  );
94 
95  const dimensionedScalar Xsolvent0
96  (
97  Ysolvent0*Wsolvent/((1 - Ysolvent0)*Wwax + Ysolvent0*Wsolvent)
98  );
99 
100  mu_ = pow(muWax_/muSolvent_, (1 - Xsolvent)/(1 - Xsolvent0))*muSolvent_;
102 }
103 
104 
105 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
106 
107 waxSolventViscosity::waxSolventViscosity
108 (
110  const dictionary& dict,
111  volScalarField& mu
112 )
113 :
114  filmViscosityModel(typeName, film, dict, mu),
115  muWax_
116  (
117  IOobject
118  (
119  typeName + ":muWax",
120  film.regionMesh().time().timeName(),
121  film.regionMesh(),
124  ),
125  film.regionMesh(),
127  zeroGradientFvPatchScalarField::typeName
128  ),
130  (
132  (
133  film,
134  coeffDict_.subDict("muWax"),
135  muWax_
136  )
137  ),
138  muSolvent_
139  (
140  IOobject
141  (
142  typeName + ":muSolvent",
143  film.regionMesh().time().timeName(),
144  film.regionMesh(),
147  ),
148  film.regionMesh(),
150  zeroGradientFvPatchScalarField::typeName
151  ),
153  (
155  (
156  film,
157  coeffDict_.subDict("muSolvent"),
158  muSolvent_
159  )
160  )
161 {}
162 
163 
164 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
165 
167 {}
168 
169 
170 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
171 
173 (
174  const volScalarField& p,
175  const volScalarField& T
176 )
177 {
178  muWaxModel_->correct(p, T);
179  muSolventModel_->correct(p, T);
180 
181  correctMu();
182 }
183 
184 
185 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
186 
187 } // End namespace surfaceFilmModels
188 } // End namespace regionModels
189 } // End namespace Foam
190 
191 // ************************************************************************* //
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
static autoPtr< filmViscosityModel > New(surfaceFilmRegionModel &film, const dictionary &dict, volScalarField &mu)
Return a reference to the selected phase change model.
UniformDimensionedField< scalar > uniformDimensionedScalarField
const dimensionSet dimDynamicViscosity
const dictionary coeffDict_
Coefficients dictionary.
Definition: subModelBase.H:79
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
Macros for easy insertion into run-time selection tables.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:692
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
autoPtr< filmViscosityModel > muWaxModel_
Wax viscosity model.
Base class for surface film viscosity models.
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:61
autoPtr< filmViscosityModel > muSolventModel_
Solvent viscosity model.
volScalarField & mu_
Reference to the viscosity field.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual void correct(const volScalarField &p, const volScalarField &T)
Correct.
void correctBoundaryConditions()
Correct boundary field.
volScalarField & p
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
defineTypeNameAndDebug(kinematicSingleLayer, 0)
Namespace for OpenFOAM.