Lee.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) 2017-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::meltingEvaporationModels::Lee
28
29Description
30 Mass transfer Lee model. Simple model driven by field value difference as:
31
32 \f[
33 \dot{m} = C \rho \alpha (T - T_{activate})/T_{activate}
34 \f]
35
36 where C is a model constant.
37
38 if C > 0:
39 \f[
40 \dot{m} = C \rho \alpha (T - T_{activate})/T_{activate}
41 \f]
42 for \f[ T > T_{activate} \f]
43
44 and
45
46 \f[ mDot = 0.0 \f] for \f[ T < T_{activate} \f]
47
48
49 if C < 0:
50 \f[
51 \dot{m} = -C \rho \alpha (T_{activate} - T)/T_{activate}
52 \f]
53 for \f[ T < T_{activate} \f]
54
55 and
56 \f[ \dot{m} = 0.0 \f] for \f[ T > T_{activate} \f]
57
58 Based on the reference:
59 -# W. H. Lee. "A Pressure Iteration Scheme for Two-Phase Modeling".
60 Technical Report LA-UR 79-975. Los Alamos Scientific Laboratory,
61 Los Alamos, New Mexico. 1979.
62
63Usage
64 Example usage:
65 \verbatim
66 massTransferModel
67 (
68 (solid to liquid)
69 {
70 type Lee;
71 C 40;
72 Tactivate 302.78;
73 }
74 );
75 \endverbatim
76
77 Where:
78
79 \table
80 Property | Description | Required | Default value
81 Tactivate | Activation temperature | yes
82 C | Model constant | yes
83 includeVolChange | Volumen change | no | yes
84 species | Specie name on the other phase | no | none
85 \endtable
86
87SourceFiles
88 Lee.C
89
90\*---------------------------------------------------------------------------*/
91
92#ifndef meltingEvaporationModels_Lee_H
93#define meltingEvaporationModels_Lee_H
94
95#include "InterfaceCompositionModel.H"
96
97// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *//
98
99namespace Foam
100{
101namespace meltingEvaporationModels
102{
103
104/*---------------------------------------------------------------------------*\
105 Class Lee Declaration
106\*---------------------------------------------------------------------------*/
107
108template<class Thermo, class OtherThermo>
109class Lee
110:
111 public InterfaceCompositionModel<Thermo, OtherThermo>
112{
113 // Private Data
114
115 //- Condensation coefficient [1/s]
117
118 //- Phase transition temperature
119 const dimensionedScalar Tactivate_;
120
121 //- Phase minimum value for activation
122 scalar alphaMin_;
123
124
125public:
126
127 //- Runtime type information
128 TypeName("Lee");
129
130
131 // Constructors
132
133 //- Construct from components
134 Lee
135 (
136 const dictionary& dict,
137 const phasePair& pair
138 );
139
140
141 //- Destructor
142 virtual ~Lee() = default;
143
144
145 // Member Functions
146
147 //- Explicit total mass transfer coefficient
149 (
151 );
152
153 //- Implicit mass transfer coefficient
155 (
156 label modelVariable,
157 const volScalarField& field
158 );
159
160 //- Explicit mass transfer coefficient
162 (
163 label modelVariable,
165 );
166
167 //- Return T transition between phases
168 virtual const dimensionedScalar& Tactivate() const noexcept
169 {
170 return Tactivate_;
171 }
172
173 //- Add/subtract alpha*div(U) as a source term
174 //- for alpha, substituting div(U) = mDot(1/rho1 - 1/rho2)
175 virtual bool includeDivU() const noexcept
176 {
177 return true;
178 }
179};
180
181
182// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
183
184} // End namespace meltingEvaporationModels
185} // End namespace Foam
186
187// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188
189#ifdef NoRepository
190# include "Lee.C"
191#endif
192
193// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
194
195#endif
196
197// ************************************************************************* //
Base class for interface composition models, templated on the two thermodynamic models either side of...
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mass transfer Lee model. Simple model driven by field value difference as:
Definition: Lee.H:134
TypeName("Lee")
Runtime type information.
virtual ~Lee()=default
Destructor.
virtual tmp< volScalarField > Kexp(const volScalarField &field)
Explicit total mass transfer coefficient.
Definition: Lee.C:53
virtual tmp< volScalarField > KSu(label modelVariable, const volScalarField &field)
Explicit mass transfer coefficient.
Definition: Lee.C:134
virtual const dimensionedScalar & Tactivate() const noexcept
Return T transition between phases.
Definition: Lee.H:190
virtual bool includeDivU() const noexcept
Definition: Lee.H:197
virtual tmp< volScalarField > KSp(label modelVariable, const volScalarField &field)
Implicit mass transfer coefficient.
Definition: Lee.C:91
const phasePair & pair() const
The phase pair.
modelVariable
Enumeration for variable based mass transfer models.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:56
A class for managing temporary objects.
Definition: tmp.H:65
rDeltaTY field()
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const direction noexcept
Definition: Scalar.H:223
dictionary dict
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73