linearTsubDiameter.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) 2018-2023 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 "linearTsubDiameter.H"
27 #include "phaseSystem.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace diameterModels
36 {
39 }
40 }
41 
42 
43 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 
46 (
47  const dictionary& diameterProperties,
48  const phaseModel& phase
49 )
50 :
51  spherical(diameterProperties, phase),
52  liquidPhaseName_(diameterProperties.lookup("liquidPhase")),
53  d2_("d2", dimLength, diameterProperties.lookupOrDefault("d2", 0.0015)),
54  Tsub2_
55  (
56  "Tsub2",
58  diameterProperties.lookupOrDefault("Tsub2", 0)
59  ),
60  d1_("d1", dimLength, diameterProperties.lookupOrDefault("d1", 0.00015)),
61  Tsub1_
62  (
63  "Tsub1",
65  diameterProperties.lookupOrDefault("Tsub1", 13.5)
66  ),
67  d_
68  (
69  IOobject
70  (
71  IOobject::groupName("d", phase.name()),
72  phase.time().name(),
73  phase.mesh()
74  ),
75  phase.mesh(),
76  d1_
77  )
78 {
79  Info<< " d2: " << d2_.value() << endl
80  << " Tsub2: " << Tsub2_.value() << endl
81  << " d1: " << d1_.value() << endl
82  << " Tsub1: " << Tsub1_.value() << endl;
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87 
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
93 
95 {
96  return d_;
97 }
98 
99 
101 {
102  const phaseSystem& fluid = phase().fluid();
103  const phaseModel& liquid = fluid.phases()[liquidPhaseName_];
104  const phaseInterface interface(phase(), liquid);
105 
106  if
107  (
109  <
111  >(interface)
112  )
113  {
114  const interfaceSaturationTemperatureModel& satModel =
116  <
118  >(interface);
119 
120  const volScalarField Tsub
121  (
122  satModel.Tsat(liquid.fluidThermo().p()) - liquid.thermo().T()
123  );
124 
125  d_ = max
126  (
127  d1_,
128  min
129  (
130  d2_,
131  (d1_*(Tsub - Tsub2_) + d2_*(Tsub - Tsub1_))/(Tsub2_ - Tsub1_)
132  )
133  );
134  }
135 }
136 
137 
139 {
141 
142  diameterProperties().lookup("liquidPhase") >> liquidPhaseName_;
143 
144  d2_ =
146  (
147  dimLength,
148  diameterProperties().lookupOrDefault<scalar>("d2", 0.0015)
149  );
150  Tsub2_ =
152  (
154  diameterProperties().lookupOrDefault<scalar>("Tsub2", 0)
155  );
156  d1_ =
158  (
159  dimLength,
160  diameterProperties().lookupOrDefault("d1", 0.00015)
161  );
162  Tsub1_ =
164  (
166  diameterProperties().lookupOrDefault<scalar>("Tsub1", 13.5)
167  );
168 
169  return true;
170 }
171 
172 
173 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
Abstract base-class for dispersed-phase particle diameter models.
Definition: diameterModel.H:52
virtual bool read(const dictionary &phaseProperties)
Read phaseProperties dictionary.
Definition: diameterModel.C:62
linearTsub(const dictionary &diameterProperties, const phaseModel &phase)
Construct from dictionary and phase.
virtual void correct()
Correct the model.
virtual tmp< volScalarField > d() const
Get the diameter field.
virtual bool read(const dictionary &phaseProperties)
Read phaseProperties dictionary.
Base class for models which represent spherical diameter models, providing a common implementation of...
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const Type & value() const
Return const reference to value.
Wrapper around saturationTemperatureModel to facilitate convenient construction on interfaces.
tmp< volScalarField > Tsat(const volScalarField &p) const
Saturation temperature.
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:53
Class to represent an interface between phases. Derivations can further specify the configuration of ...
Helper class to manage multi-specie phase properties.
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
bool foundInterfacialModel(const phaseInterface &interface) const
Check availability of a sub model for a given interface.
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:102
const ModelType & lookupInterfacialModel(const phaseInterface &interface) const
Return a sub model for an interface.
A class for managing temporary objects.
Definition: tmp.H:55
addToRunTimeSelectionTable(diameterModel, constant, dictionary)
defineTypeNameAndDebug(constant, 0)
Namespace for OpenFOAM.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
const dimensionSet dimLength
const dimensionSet dimTemperature
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.