standardPhaseChange.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-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 "standardPhaseChange.H"
28 #include "thermoSingleLayer.H"
29 #include "zeroField.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace regionModels
36 {
37 namespace surfaceFilmModels
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 defineTypeNameAndDebug(standardPhaseChange, 0);
43 
45 (
46  phaseChangeModel,
47  standardPhaseChange,
48  dictionary
49 );
50 
51 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
52 
54 (
55  const scalar Re,
56  const scalar Sc
57 ) const
58 {
59  if (Re < 5.0E+05)
60  {
61  return 0.664*sqrt(Re)*cbrt(Sc);
62  }
63  else
64  {
65  return 0.037*pow(Re, 0.8)*cbrt(Sc);
66  }
67 }
68 
69 
70 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
71 
72 standardPhaseChange::standardPhaseChange
73 (
75  const dictionary& dict
76 )
77 :
78  phaseChangeModel(typeName, film, dict),
79  deltaMin_(readScalar(coeffDict_.lookup("deltaMin"))),
80  L_(readScalar(coeffDict_.lookup("L"))),
81  TbFactor_(coeffDict_.lookupOrDefault<scalar>("TbFactor", 1.1)),
82  YInfZero_(coeffDict_.lookupOrDefault<Switch>("YInfZero", false))
83 {}
84 
85 
86 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87 
89 {}
90 
91 
92 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
93 
94 template<class YInfType>
96 (
97  const scalar dt,
98  scalarField& availableMass,
99  scalarField& dMass,
100  scalarField& dEnergy,
101  YInfType YInf
102 )
103 {
104  const thermoSingleLayer& film = filmType<thermoSingleLayer>();
105 
106  // Set local thermo properties
107  const SLGThermo& thermo = film.thermo();
108  const filmThermoModel& filmThermo = film.filmThermo();
109  const label vapId = thermo.carrierId(filmThermo.name());
110 
111  // Retrieve fields from film model
112  const scalarField& delta = film.delta();
113  const scalarField& pInf = film.pPrimary();
114  const scalarField& T = film.T();
115  const scalarField& hs = film.hs();
116  const scalarField& rho = film.rho();
117  const scalarField& rhoInf = film.rhoPrimary();
118  const scalarField& muInf = film.muPrimary();
119  const scalarField& magSf = film.magSf();
120  const vectorField dU(film.UPrimary() - film.Us());
121  const scalarField limMass
122  (
123  max(scalar(0), availableMass - deltaMin_*rho*magSf)
124  );
125 
126  // Molecular weight of vapour [kg/kmol]
127  const scalar Wvap = thermo.carrier().W(vapId);
128 
129  // Molecular weight of liquid [kg/kmol]
130  const scalar Wliq = filmThermo.W();
131 
132  forAll(dMass, celli)
133  {
134  scalar dm = 0;
135 
136  if (delta[celli] > deltaMin_)
137  {
138  // Cell pressure [Pa]
139  const scalar pc = pInf[celli];
140 
141  // Calculate the boiling temperature
142  const scalar Tb = filmThermo.Tb(pc);
143 
144  // Local temperature - impose lower limit of 200 K for stability
145  const scalar Tloc = min(TbFactor_*Tb, max(200.0, T[celli]));
146 
147  // Saturation pressure [Pa]
148  const scalar pSat = filmThermo.pv(pc, Tloc);
149 
150  // Latent heat [J/kg]
151  const scalar hVap = filmThermo.hl(pc, Tloc);
152 
153  // Calculate mass transfer
154  if (pSat >= 0.95*pc)
155  {
156  // Boiling
157  const scalar Cp = filmThermo.Cp(pc, Tloc);
158  const scalar Tcorr = max(0.0, T[celli] - Tb);
159  const scalar qCorr = limMass[celli]*Cp*(Tcorr);
160  dm = qCorr/hVap;
161  }
162  else
163  {
164  // Primary region density [kg/m3]
165  const scalar rhoInfc = rhoInf[celli];
166 
167  // Primary region viscosity [Pa.s]
168  const scalar muInfc = muInf[celli];
169 
170  // Reynolds number
171  const scalar Re = rhoInfc*mag(dU[celli])*L_/muInfc;
172 
173  // Vapour mass fraction at interface
174  const scalar Ys = Wliq*pSat/(Wliq*pSat + Wvap*(pc - pSat));
175 
176  // Vapour diffusivity [m2/s]
177  const scalar Dab = filmThermo.D(pc, Tloc);
178 
179  // Schmidt number
180  const scalar Sc = muInfc/(rhoInfc*(Dab + rootVSmall));
181 
182  // Sherwood number
183  const scalar Sh = this->Sh(Re, Sc);
184 
185  // Mass transfer coefficient [m/s]
186  const scalar hm = Sh*Dab/(L_ + rootVSmall);
187 
188  // Add mass contribution to source
189  dm = dt*magSf[celli]*rhoInfc*hm*(Ys - YInf[celli])/(1.0 - Ys);
190  }
191 
192  dMass[celli] += min(limMass[celli], max(dm, 0));
193 
194  // Heat is assumed to be removed by heat-transfer to the wall
195  // so the energy remains unchanged by the phase-change.
196  dEnergy[celli] += dm*hs[celli];
197  // dEnergy[celli] += dm*(hs[celli] + hVap);
198  }
199  }
200 }
201 
202 
204 (
205  const scalar dt,
206  scalarField& availableMass,
207  scalarField& dMass,
208  scalarField& dEnergy
209 )
210 {
211  if (YInfZero_)
212  {
213  correctModel(dt, availableMass, dMass, dEnergy, zeroField());
214  }
215  else
216  {
217  const thermoSingleLayer& film = filmType<thermoSingleLayer>();
218  const label vapId = film.thermo().carrierId(film.filmThermo().name());
219  const scalarField& YInf = film.YPrimary()[vapId];
220 
221  correctModel(dt, availableMass, dMass, dEnergy, YInf);
222  }
223 }
224 
225 
226 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
227 
228 } // End namespace surfaceFilmModels
229 } // End namespace regionModels
230 } // End namespace Foam
231 
232 // ************************************************************************* //
scalar delta
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
virtual const word & name() const =0
Return the specie name.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual scalar Cp(const scalar p, const scalar T) const =0
Return specific heat capacity [J/kg/K].
virtual scalar hl(const scalar p, const scalar T) const =0
Return latent heat [J/kg].
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual const volVectorField & Us() const
Return the film surface velocity [m/s].
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
const volScalarField & rhoPrimary() const
Density [kg/m3].
rhoReactionThermo & thermo
Definition: createFields.H:28
addToRunTimeSelectionTable(surfaceFilmRegionModel, kinematicSingleLayer, mesh)
Macros for easy insertion into run-time selection tables.
Base class for surface film phase change models.
const volScalarField & muPrimary() const
Viscosity [Pa.s].
const filmThermoModel & filmThermo() const
Film thermo.
const volVectorField & UPrimary() const
Velocity [m/s].
scalar Sh(const scalar Re, const scalar Sc) const
Return Sherwood number as a function of Reynolds and Schmidt numbers.
A class representing the concept of a field of 0 used to avoid unnecessary manipulations for objects ...
Definition: zeroField.H:50
const scalar TbFactor_
Boiling temperature factor / [].
dimensionedScalar cbrt(const dimensionedScalar &ds)
const PtrList< volScalarField > & YPrimary() const
Specie mass fractions [0-1].
Thermo package for (S)olids (L)iquids and (G)ases Takes reference to thermo package, and provides:
Definition: SLGThermo.H:62
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
virtual scalar D(const scalar p, const scalar T) const =0
Return diffusivity [m2/s].
virtual scalar Tb(const scalar p) const =0
Return boiling temperature [K].
const volScalarField & pPrimary() const
Pressure [Pa].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & delta() const
Return const access to the film thickness [m].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void correctModel(const scalar dt, scalarField &availableMass, scalarField &dMass, scalarField &dEnergy, YInfType YInf)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual scalar W() const =0
Return molecular weight [kg/kmol].
label carrierId(const word &cmptName, bool allowNotFound=false) const
Index of carrier component.
Definition: SLGThermo.C:148
virtual const volScalarField & magSf() const
Return the face area magnitudes / [m2].
virtual const volScalarField & hs() const
Return the film sensible enthalpy [J/kg].
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual const volScalarField & rho() const
Return the film density [kg/m3].
const scalar deltaMin_
Minimum film height for model to be active.
defineTypeNameAndDebug(kinematicSingleLayer, 0)
const SLGThermo & thermo() const
Return const reference to the SLGThermo object.
virtual scalar pv(const scalar p, const scalar T) const =0
Return vapour pressure [Pa].
Thermodynamic form of single-cell layer surface film model.
Namespace for OpenFOAM.
virtual const volScalarField & T() const
Return the film mean temperature [K].