HeatAndMassTransferPhaseSystem.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) 2015-2016 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 
27 
28 #include "BlendedInterfacialModel.H"
29 #include "heatTransferModel.H"
30 #include "massTransferModel.H"
31 
32 #include "HashPtrTable.H"
33 
34 #include "fvcDiv.H"
35 #include "fvmSup.H"
36 #include "fvMatrix.H"
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
41 template<class BasePhaseSystem>
44 (
45  const fvMesh& mesh
46 )
47 :
48  BasePhaseSystem(mesh)
49 {
51  (
52  "heatTransfer",
53  heatTransferModels_
54  );
55 
57  (
58  "massTransfer",
59  massTransferModels_
60  );
61 
63  (
65  this->phasePairs_,
66  phasePairIter
67  )
68  {
69  const phasePair& pair(phasePairIter());
70 
71  if (pair.ordered())
72  {
73  continue;
74  }
75 
76  // Initialy assume no mass transfer
77 
78  dmdt_.insert
79  (
80  pair,
81  new volScalarField
82  (
83  IOobject
84  (
85  IOobject::groupName("dmdt", pair.name()),
86  this->mesh().time().timeName(),
87  this->mesh(),
90  ),
91  this->mesh(),
93  )
94  );
95 
96  dmdtExplicit_.insert
97  (
98  pair,
99  new volScalarField
100  (
101  IOobject
102  (
103  IOobject::groupName("dmdtExplicit", pair.name()),
104  this->mesh().time().timeName(),
105  this->mesh()
106  ),
107  this->mesh(),
109  )
110  );
111 
112  volScalarField H1(heatTransferModels_[pair][pair.first()]->K());
113  volScalarField H2(heatTransferModels_[pair][pair.second()]->K());
114 
115  Tf_.insert
116  (
117  pair,
118  new volScalarField
119  (
120  IOobject
121  (
122  IOobject::groupName("Tf", pair.name()),
123  this->mesh().time().timeName(),
124  this->mesh(),
127  ),
128  (
129  H1*pair.phase1().thermo().T()
130  + H2*pair.phase2().thermo().T()
131  )
132  /max
133  (
134  H1 + H2,
136  ),
137  zeroGradientFvPatchScalarField::typeName
138  )
139  );
140  Tf_[pair]->correctBoundaryConditions();
141  }
142 }
143 
144 
145 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
146 
147 template<class BasePhaseSystem>
150 {}
151 
152 
153 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
154 
155 template<class BasePhaseSystem>
157 (
158  const phaseModel& phase
159 ) const
160 {
161  return true;
162 }
163 
164 
165 template<class BasePhaseSystem>
168 (
169  const phasePairKey& key
170 ) const
171 {
172  const scalar dmdtSign(Pair<word>::compare(dmdt_.find(key).key(), key));
173 
174  return dmdtSign**dmdt_[key];
175 }
176 
177 
178 template<class BasePhaseSystem>
181 (
182  const Foam::phaseModel& phase
183 ) const
184 {
185  tmp<volScalarField> tdmdt
186  (
187  new volScalarField
188  (
189  IOobject
190  (
191  IOobject::groupName("dmdt", phase.name()),
192  this->mesh_.time().timeName(),
193  this->mesh_
194  ),
195  this->mesh_,
197  )
198  );
199 
201  (
203  this->phasePairs_,
204  phasePairIter
205  )
206  {
207  const phasePair& pair(phasePairIter());
208 
209  if (pair.ordered())
210  {
211  continue;
212  }
213 
214  const phaseModel* phase1 = &pair.phase1();
215  const phaseModel* phase2 = &pair.phase2();
216 
217  forAllConstIter(phasePair, pair, iter)
218  {
219  if (phase1 == &phase)
220  {
221  tdmdt.ref() += this->dmdt(pair);
222  }
223 
224  Swap(phase1, phase2);
225  }
226  }
227 
228  return tdmdt;
229 }
230 
231 
232 template<class BasePhaseSystem>
235 {
236  autoPtr<phaseSystem::momentumTransferTable>
238 
239  phaseSystem::momentumTransferTable& eqns = eqnsPtr();
240 
241  // Source term due to mass trasfer
243  (
245  this->phasePairs_,
246  phasePairIter
247  )
248  {
249  const phasePair& pair(phasePairIter());
250 
251  if (pair.ordered())
252  {
253  continue;
254  }
255 
256  const volVectorField& U1(pair.phase1().U());
257  const volVectorField& U2(pair.phase2().U());
258 
259  const volScalarField dmdt(this->dmdt(pair));
260  const volScalarField dmdt21(posPart(dmdt));
261  const volScalarField dmdt12(negPart(dmdt));
262 
263  *eqns[pair.phase1().name()] += dmdt21*U2 - fvm::Sp(dmdt21, U1);
264  *eqns[pair.phase2().name()] -= dmdt12*U1 - fvm::Sp(dmdt12, U2);
265  }
266 
267  return eqnsPtr;
268 }
269 
270 
271 template<class BasePhaseSystem>
274 {
275  autoPtr<phaseSystem::heatTransferTable> eqnsPtr
276  (
278  );
279 
280  phaseSystem::heatTransferTable& eqns = eqnsPtr();
281 
282  forAll(this->phaseModels_, phasei)
283  {
284  const phaseModel& phase = this->phaseModels_[phasei];
285 
286  eqns.insert
287  (
288  phase.name(),
289  new fvScalarMatrix(phase.thermo().he(), dimEnergy/dimTime)
290  );
291  }
292 
293  // Heat transfer with the interface
295  (
296  heatTransferModelTable,
297  heatTransferModels_,
298  heatTransferModelIter
299  )
300  {
301  const phasePair& pair
302  (
303  this->phasePairs_[heatTransferModelIter.key()]
304  );
305 
306  const phaseModel* phase = &pair.phase1();
307  const phaseModel* otherPhase = &pair.phase2();
308 
309  const volScalarField& Tf(*Tf_[pair]);
310 
311  const volScalarField K1
312  (
313  heatTransferModelIter()[pair.first()]->K()
314  );
315  const volScalarField K2
316  (
317  heatTransferModelIter()[pair.second()]->K()
318  );
319  const volScalarField KEff
320  (
321  K1*K2
322  /max
323  (
324  K1 + K2,
326  )
327  );
328 
329  const volScalarField* K = &K1;
330  const volScalarField* otherK = &K2;
331 
332  forAllConstIter(phasePair, pair, iter)
333  {
334  const volScalarField& he(phase->thermo().he());
335  volScalarField Cpv(phase->thermo().Cpv());
336 
337  *eqns[phase->name()] +=
338  (*K)*(Tf - phase->thermo().T())
339  + KEff/Cpv*he - fvm::Sp(KEff/Cpv, he);
340 
341  Swap(phase, otherPhase);
342  Swap(K, otherK);
343  }
344  }
345 
346  // Source term due to mass transfer
348  (
350  this->phasePairs_,
351  phasePairIter
352  )
353  {
354  const phasePair& pair(phasePairIter());
355 
356  if (pair.ordered())
357  {
358  continue;
359  }
360 
361  const phaseModel& phase1 = pair.phase1();
362  const phaseModel& phase2 = pair.phase2();
363 
364  const volScalarField& he1(phase1.thermo().he());
365  const volScalarField& he2(phase2.thermo().he());
366 
367  const volScalarField& K1(phase1.K());
368  const volScalarField& K2(phase2.K());
369 
370  const volScalarField dmdt(this->dmdt(pair));
371  const volScalarField dmdt21(posPart(dmdt));
372  const volScalarField dmdt12(negPart(dmdt));
373  const volScalarField& Tf(*Tf_[pair]);
374 
375  *eqns[phase1.name()] +=
376  dmdt21*(phase1.thermo().he(phase1.thermo().p(), Tf))
377  - fvm::Sp(dmdt21, he1)
378  + dmdt21*(K2 - K1);
379 
380  *eqns[phase2.name()] -=
381  dmdt12*(phase2.thermo().he(phase2.thermo().p(), Tf))
382  - fvm::Sp(dmdt12, he2)
383  + dmdt12*(K1 - K2);
384  }
385 
386  return eqnsPtr;
387 }
388 
389 
390 template<class BasePhaseSystem>
392 {
393  if (BasePhaseSystem::read())
394  {
395  bool readOK = true;
396 
397  // Models ...
398 
399  return readOK;
400  }
401  else
402  {
403  return false;
404  }
405 }
406 
407 
408 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:56
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
label phasei
Definition: pEqn.H:27
#define K1
Definition: SHA1.C:167
HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash > phasePairTable
Definition: phaseSystem.H:133
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
virtual autoPtr< phaseSystem::momentumTransferTable > momentumTransfer() const
Return the momentum transfer matrices.
#define K2
Definition: SHA1.C:168
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
phasePairTable phasePairs_
Phase pairs.
Definition: phaseSystem.H:167
static int compare(const Pair< word > &a, const Pair< word > &b)
Compare Pairs.
Definition: Pair.H:141
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:114
virtual bool transfersMass(const phaseModel &phase) const
Return true if there is mass transfer for phase.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
dimensionedScalar posPart(const dimensionedScalar &ds)
phaseModelList phaseModels_
Phase models.
Definition: phaseSystem.H:164
void insert(const word &, T *)
Add at head of dictionary.
static const dimensionSet dimK
Coefficient dimensions.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
HashPtrTable< fvScalarMatrix, word, string::hash > heatTransferTable
Definition: phaseSystem.H:109
bool read(const char *, int32_t &)
Definition: int32IO.C:85
void Swap(T &a, T &b)
Definition: Swap.H:43
const fvMesh & mesh() const
Constant access the mesh.
Definition: phaseSystemI.H:28
phaseSystem::momentumTransferTable & momentumTransfer(momentumTransferPtr())
static word groupName(Name name, const word &group)
virtual ~HeatAndMassTransferPhaseSystem()
Destructor.
HashPtrTable< fvVectorMatrix, word, string::hash > momentumTransferTable
Definition: phaseSystem.H:100
const fvMesh & mesh_
Reference to the mesh.
Definition: phaseSystem.H:161
volVectorField & U1
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the interfacial mass flow rate.
Calculate the divergence of the given field.
volScalarField & he2
Definition: EEqns.H:3
HeatAndMassTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
const dimensionSet dimEnergy
const dimensionSet dimDensity
const word & name() const
Definition: phaseModel.H:109
void generatePairsAndSubModels(const word &modelName, HashTable< autoPtr< modelType >, phasePairKey, phasePairKey::hash > &models)
Generate pairs and sub-model tables.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
volVectorField & U2
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
void correctBoundaryConditions()
Correct boundary field.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
A class for managing temporary objects.
Definition: PtrList.H:54
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
virtual bool read()
Read base phaseProperties dictionary.
Calculate the matrix for implicit and explicit sources.
dimensionedScalar negPart(const dimensionedScalar &ds)
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243