Lee.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2017-2022 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
26\*---------------------------------------------------------------------------*/
27
28#include "Lee.H"
30
31
32// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33
34template<class Thermo, class OtherThermo>
36(
37 const dictionary& dict,
38 const phasePair& pair
39)
40:
41 InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
42 C_("C", inv(dimTime), dict),
43 Tactivate_("Tactivate", dimTemperature, dict),
44 alphaMin_(dict.getOrDefault<scalar>("alphaMin", 0))
45{}
46
47
48// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
49
50template<class Thermo, class OtherThermo>
53(
54 const volScalarField& refValue
55)
56{
57 {
58 const volScalarField from
59 (
60 min(max(this->pair().from(), scalar(0)), scalar(1))
61 );
62
63 const volScalarField coeff
64 (
65 C_*from*this->pair().from().rho()*pos(from - alphaMin_)
66 *(refValue - Tactivate_)
67 /Tactivate_
68 );
69
70 if (sign(C_.value()) > 0)
71 {
72 return
73 (
74 coeff*pos(refValue - Tactivate_)
75 );
76 }
77 else
78 {
79 return
80 (
81 coeff*pos(Tactivate_ - refValue)
82 );
83 }
84 }
85}
86
87
88template<class Thermo, class OtherThermo>
91(
92 label variable,
93 const volScalarField& refValue
94)
95{
96 if (this->modelVariable_ == variable)
97 {
99 (
100 min(max(this->pair().from(), scalar(0)), scalar(1))
101 );
102
103 const volScalarField coeff
104 (
105 C_*from*this->pair().from().rho()*pos(from - alphaMin_)
106 /Tactivate_
107 );
108
109 if (sign(C_.value()) > 0)
110 {
111 return
112 (
113 coeff*pos(refValue - Tactivate_)
114 );
115 }
116 else
117 {
118 return
119 (
120 -coeff*pos(Tactivate_ - refValue)
121 );
122 }
123 }
124 else
125 {
126 return tmp<volScalarField> ();
127 }
128}
129
130
131template<class Thermo, class OtherThermo>
134(
135 label variable,
136 const volScalarField& refValue
137)
138{
139 if (this->modelVariable_ == variable)
140 {
141 volScalarField from
142 (
143 min(max(this->pair().from(), scalar(0)), scalar(1))
144 );
145
146 const volScalarField coeff
147 (
148 C_*from*this->pair().from().rho()*pos(from - alphaMin_)
149 );
150
151 if (sign(C_.value()) > 0)
152 {
153 return
154 (
155 -coeff*pos(refValue - Tactivate_)
156 );
157 }
158 else
159 {
160 return
161 (
162 -coeff*pos(Tactivate_ - refValue)
163 );
164 }
165 }
166 else
167 {
168 return tmp<volScalarField> ();
169 }
170}
171
172
173// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
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
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 tmp< volScalarField > KSp(label modelVariable, const volScalarField &field)
Implicit mass transfer coefficient.
Definition: Lee.C:91
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
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar sign(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:54
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
dictionary dict