LiquidEvaporation.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 | Copyright (C) 2011-2015 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 "LiquidEvaporation.H"
27 #include "specie.H"
28 #include "mathematicalConstants.H"
29 
30 using namespace Foam::constant::mathematical;
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 template<class CloudType>
36 (
37  const label cellI
38 ) const
39 {
40  scalarField Xc(this->owner().thermo().carrier().Y().size());
41 
42  forAll(Xc, i)
43  {
44  Xc[i] =
45  this->owner().thermo().carrier().Y()[i][cellI]
46  /this->owner().thermo().carrier().W(i);
47  }
48 
49  return Xc/sum(Xc);
50 }
51 
52 
53 template<class CloudType>
55 (
56  const scalar Re,
57  const scalar Sc
58 ) const
59 {
60  return 2.0 + 0.6*Foam::sqrt(Re)*cbrt(Sc);
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
66 template<class CloudType>
68 (
69  const dictionary& dict,
70  CloudType& owner
71 )
72 :
73  PhaseChangeModel<CloudType>(dict, owner, typeName),
74  liquids_(owner.thermo().liquids()),
75  activeLiquids_(this->coeffDict().lookup("activeLiquids")),
76  liqToCarrierMap_(activeLiquids_.size(), -1),
77  liqToLiqMap_(activeLiquids_.size(), -1)
78 {
79  if (activeLiquids_.size() == 0)
80  {
81  WarningIn
82  (
83  "Foam::LiquidEvaporation<CloudType>::LiquidEvaporation"
84  "("
85  "const dictionary& dict, "
86  "CloudType& owner"
87  ")"
88  ) << "Evaporation model selected, but no active liquids defined"
89  << nl << endl;
90  }
91  else
92  {
93  Info<< "Participating liquid species:" << endl;
94 
95  // Determine mapping between liquid and carrier phase species
96  forAll(activeLiquids_, i)
97  {
98  Info<< " " << activeLiquids_[i] << endl;
99  liqToCarrierMap_[i] =
100  owner.composition().carrierId(activeLiquids_[i]);
101  }
102 
103  // Determine mapping between model active liquids and global liquids
104  const label idLiquid = owner.composition().idLiquid();
105  forAll(activeLiquids_, i)
106  {
107  liqToLiqMap_[i] =
108  owner.composition().localId(idLiquid, activeLiquids_[i]);
109  }
110  }
111 }
112 
113 
114 template<class CloudType>
116 (
118 )
119 :
121  liquids_(pcm.owner().thermo().liquids()),
122  activeLiquids_(pcm.activeLiquids_),
123  liqToCarrierMap_(pcm.liqToCarrierMap_),
124  liqToLiqMap_(pcm.liqToLiqMap_)
125 {}
126 
127 
128 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
129 
130 template<class CloudType>
132 {}
133 
134 
135 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
136 
137 template<class CloudType>
139 (
140  const scalar dt,
141  const label cellI,
142  const scalar Re,
143  const scalar Pr,
144  const scalar d,
145  const scalar nu,
146  const scalar T,
147  const scalar Ts,
148  const scalar pc,
149  const scalar Tc,
150  const scalarField& X,
151  scalarField& dMassPC
152 ) const
153 {
154  // immediately evaporate mass that has reached critical condition
155  if ((liquids_.Tc(X) - T) < SMALL)
156  {
157  if (debug)
158  {
159  WarningIn
160  (
161  "void Foam::LiquidEvaporation<CloudType>::calculate"
162  "("
163  "const scalar, "
164  "const label, "
165  "const scalar, "
166  "const scalar, "
167  "const scalar, "
168  "const scalar, "
169  "const scalar, "
170  "const scalar, "
171  "const scalar, "
172  "const scalar, "
173  "const scalarField&, "
174  "scalarField&"
175  ") const"
176  ) << "Parcel reached critical conditions: "
177  << "evaporating all avaliable mass" << endl;
178  }
179 
180  forAll(activeLiquids_, i)
181  {
182  const label lid = liqToLiqMap_[i];
183  dMassPC[lid] = GREAT;
184  }
185 
186  return;
187  }
188 
189  // construct carrier phase species volume fractions for cell, cellI
190  const scalarField Xc(calcXc(cellI));
191 
192  // calculate mass transfer of each specie in liquid
193  forAll(activeLiquids_, i)
194  {
195  const label gid = liqToCarrierMap_[i];
196  const label lid = liqToLiqMap_[i];
197 
198  // vapour diffusivity [m2/s]
199  const scalar Dab = liquids_.properties()[lid].D(pc, Ts);
200 
201  // saturation pressure for species i [pa]
202  // - carrier phase pressure assumed equal to the liquid vapour pressure
203  // close to the surface
204  // NOTE: if pSat > pc then particle is superheated
205  // calculated evaporation rate will be greater than that of a particle
206  // at boiling point, but this is not a boiling model
207  const scalar pSat = liquids_.properties()[lid].pv(pc, T);
208 
209  // Schmidt number
210  const scalar Sc = nu/(Dab + ROOTVSMALL);
211 
212  // Sherwood number
213  const scalar Sh = this->Sh(Re, Sc);
214 
215  // mass transfer coefficient [m/s]
216  const scalar kc = Sh*Dab/(d + ROOTVSMALL);
217 
218  // vapour concentration at surface [kmol/m3] at film temperature
219  const scalar Cs = pSat/(RR*Ts);
220 
221  // vapour concentration in bulk gas [kmol/m3] at film temperature
222  const scalar Cinf = Xc[gid]*pc/(RR*Ts);
223 
224  // molar flux of vapour [kmol/m2/s]
225  const scalar Ni = max(kc*(Cs - Cinf), 0.0);
226 
227  // mass transfer [kg]
228  dMassPC[lid] += Ni*pi*sqr(d)*liquids_.properties()[lid].W()*dt;
229  }
230 }
231 
232 
233 template<class CloudType>
235 (
236  const label idc,
237  const label idl,
238  const scalar p,
239  const scalar T
240 ) const
241 {
242  scalar dh = 0;
243 
244  typedef PhaseChangeModel<CloudType> parent;
245  switch (parent::enthalpyTransfer_)
246  {
247  case (parent::etLatentHeat):
248  {
249  dh = liquids_.properties()[idl].hl(p, T);
250  break;
251  }
252  case (parent::etEnthalpyDifference):
253  {
254  scalar hc = this->owner().composition().carrier().Ha(idc, p, T);
255  scalar hp = liquids_.properties()[idl].h(p, T);
256 
257  dh = hc - hp;
258  break;
259  }
260  default:
261  {
263  (
264  "Foam::scalar Foam::LiquidEvaporation<CloudType>::dh"
265  "("
266  "const label, "
267  "const label, "
268  "const scalar, "
269  "const scalar"
270  ") const"
271  ) << "Unknown enthalpyTransfer type" << abort(FatalError);
272  }
273  }
274 
275  return dh;
276 }
277 
278 
279 template<class CloudType>
281 (
282  const scalarField& X
283 ) const
284 {
285  return liquids_.Tpt(X);
286 }
287 
288 
289 template<class CloudType>
291 (
292  const scalar p,
293  const scalarField& X
294 ) const
295 {
296  return liquids_.pvInvert(p, X);
297 }
298 
299 
300 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
const CloudType & owner() const
Return const access to the owner cloud.
scalar Sh
Definition: solveChemistry.H:2
tmp< scalarField > calcXc(const label cellI) const
Calculate the carrier phase component volume fractions at cellI.
virtual scalar dh(const label idc, const label idl, const scalar p, const scalar T) const
Return the enthalpy per unit mass.
const dimensionedScalar & pSat
Definition: createFields.H:41
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
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
dimensionedScalar cbrt(const dimensionedScalar &ds)
List< label > liqToLiqMap_
Mapping between local and global liquid species.
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual scalar Tvap(const scalarField &X) const
Return vapourisation temperature.
Templated phase change model class.
Definition: ReactingCloud.H:55
dictionary dict
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
PtrList< volScalarField > & Y
Definition: createFields.H:36
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define forAll(list, i)
Definition: UList.H:421
virtual scalar TMax(const scalar p, const scalarField &X) const
Return maximum/limiting temperature.
const dictionary & properties() const
Return const access to the properties dictionary.
Definition: subModelBase.C:134
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
virtual ~LiquidEvaporation()
Destructor.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Liquid evaporation model.
const scalar RR
Universal gas constant (default in [J/(kmol K)])
error FatalError
List< label > liqToCarrierMap_
Mapping between liquid and carrier species.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68
LiquidEvaporation(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
virtual void calculate(const scalar dt, const label cellI, const scalar Re, const scalar Pr, const scalar d, const scalar nu, const scalar T, const scalar Ts, const scalar pc, const scalar Tc, const scalarField &X, scalarField &dMassPC) const
Update model.
psiReactionThermo & thermo
Definition: createFields.H:32
dimensionedSymmTensor sqr(const dimensionedVector &dv)
scalar Sh() const
Sherwood number.
List< word > activeLiquids_
List of active liquid names.
mathematical constants.
A class for managing temporary objects.
Definition: PtrList.H:118