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