ThermalPhaseChangePhaseSystem.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015-2019 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "fvcVolumeIntegrate.H"
32#include "fvmSup.H"
33
34// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
35
36template<class BasePhaseSystem>
39(
40 const phasePairKey& key
41) const
42{
43 if (!iDmdt_.found(key))
44 {
45 return phaseSystem::dmdt(key);
46 }
47
48 const scalar dmdtSign(Pair<word>::compare(iDmdt_.find(key).key(), key));
49
50 return dmdtSign**iDmdt_[key];
51}
52
53
54template<class BasePhaseSystem>
57(
58 const phasePairKey& key
59) const
60{
61 if (!wDmdt_.found(key))
62 {
63 return phaseSystem::dmdt(key);
64 }
65
66 const scalar dmdtSign(Pair<word>::compare(wDmdt_.find(key).key(), key));
67
68 return dmdtSign**wDmdt_[key];
69}
70
71
72// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
73
74template<class BasePhaseSystem>
77(
78 const fvMesh& mesh
79)
80:
81 BasePhaseSystem(mesh),
82 volatile_(this->template getOrDefault<word>("volatile", "none")),
83 saturationModel_
84 (
85 saturationModel::New(this->subDict("saturationModel"), mesh)
86 ),
87 phaseChange_(this->lookup("phaseChange"))
88{
89
91 (
93 this->phasePairs_,
94 phasePairIter
95 )
96 {
97 const phasePair& pair(phasePairIter());
98
99 if (pair.ordered())
100 {
101 continue;
102 }
103
104 // Initially assume no mass transfer
105 iDmdt_.set
106 (
107 pair,
109 (
111 (
112 IOobject::groupName("iDmdt", pair.name()),
113 this->mesh().time().timeName(),
114 this->mesh(),
117 ),
118 this->mesh(),
120 )
121 );
122
123 // Initially assume no mass transfer
124 wDmdt_.set
125 (
126 pair,
128 (
130 (
131 IOobject::groupName("wDmdt", pair.name()),
132 this->mesh().time().timeName(),
133 this->mesh(),
136 ),
137 this->mesh(),
139 )
140 );
141
142 // Initially assume no mass transfer
144 (
145 pair,
147 (
149 (
150 IOobject::groupName("wMDotL", pair.name()),
151 this->mesh().time().timeName(),
152 this->mesh(),
155 ),
156 this->mesh(),
158 )
159 );
160 }
161}
162
163
164// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
165
166template<class BasePhaseSystem>
169{}
170
171
172// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
173
174template<class BasePhaseSystem>
177{
178 return saturationModel_();
179}
180
181
182template<class BasePhaseSystem>
185(
186 const phasePairKey& key
187) const
188{
189 return BasePhaseSystem::dmdt(key) + this->iDmdt(key) + this->wDmdt(key);
190}
191
192
193template<class BasePhaseSystem>
196{
197 PtrList<volScalarField> dmdts(BasePhaseSystem::dmdts());
198
199 forAllConstIter(iDmdtTable, iDmdt_, iDmdtIter)
200 {
201 const phasePair& pair = this->phasePairs_[iDmdtIter.key()];
202 const volScalarField& iDmdt = *iDmdtIter();
203
204 this->addField(pair.phase1(), "dmdt", iDmdt, dmdts);
205 this->addField(pair.phase2(), "dmdt", - iDmdt, dmdts);
206 }
207
208 forAllConstIter(wDmdtTable, wDmdt_, wDmdtIter)
209 {
210 const phasePair& pair = this->phasePairs_[wDmdtIter.key()];
211 const volScalarField& wDmdt = *wDmdtIter();
212
213 this->addField(pair.phase1(), "dmdt", wDmdt, dmdts);
214 this->addField(pair.phase2(), "dmdt", - wDmdt, dmdts);
215 }
216
217 return dmdts;
218}
219
220
221template<class BasePhaseSystem>
224{
226 BasePhaseSystem::heatTransfer();
227
228 phaseSystem::heatTransferTable& eqns = eqnsPtr();
229
230 // Add boundary term
232 (
234 this->phasePairs_,
235 phasePairIter
236 )
237 {
238 if (this->wMDotL_.found(phasePairIter.key()))
239 {
240 const phasePair& pair(phasePairIter());
241
242 if (pair.ordered())
243 {
244 continue;
245 }
246
247 const phaseModel& phase1 = pair.phase1();
248 const phaseModel& phase2 = pair.phase2();
249
250 *eqns[phase1.name()] += negPart(*this->wMDotL_[pair]);
251 *eqns[phase2.name()] -= posPart(*this->wMDotL_[pair]);
252
253 if
254 (
255 phase1.thermo().he().member() == "e"
256 || phase2.thermo().he().member() == "e"
257 )
258 {
259 const volScalarField dmdt
260 (
261 this->iDmdt(pair) + this->wDmdt(pair)
262 );
263
264 if (phase1.thermo().he().member() == "e")
265 {
266 *eqns[phase1.name()] +=
267 phase1.thermo().p()*dmdt/phase1.thermo().rho();
268 }
269
270 if (phase2.thermo().he().member() == "e")
271 {
272 *eqns[phase2.name()] -=
273 phase2.thermo().p()*dmdt/phase2.thermo().rho();
274 }
275 }
276 }
277 }
278
279 return eqnsPtr;
280}
281
282
283template<class BasePhaseSystem>
286{
288 BasePhaseSystem::massTransfer();
289
290 phaseSystem::massTransferTable& eqns = eqnsPtr();
291
293 (
295 this->phasePairs_,
296 phasePairIter
297 )
298 {
299 const phasePair& pair(phasePairIter());
300
301 if (pair.ordered())
302 {
303 continue;
304 }
305
306 const phaseModel& phase = pair.phase1();
307 const phaseModel& otherPhase = pair.phase2();
308
309 const PtrList<volScalarField>& Yi = phase.Y();
310
311 forAll(Yi, i)
312 {
313 if (Yi[i].member() != volatile_)
314 {
315 continue;
316 }
317
318 const word name
319 (
320 IOobject::groupName(volatile_, phase.name())
321 );
322
323 const word otherName
324 (
325 IOobject::groupName(volatile_, otherPhase.name())
326 );
327
328 // Note that the phase YiEqn does not contain a continuity error
329 // term, so these additions represent the entire mass transfer
330
331 const volScalarField dmdt(this->iDmdt(pair) + this->wDmdt(pair));
332
333 *eqns[name] += dmdt;
334 *eqns[otherName] -= dmdt;
335 }
336 }
337
338 return eqnsPtr;
339}
340
341
342template<class BasePhaseSystem>
343void
345{
347 alphatPhaseChangeWallFunction;
348
350 (
351 typename BasePhaseSystem::heatTransferModelTable,
352 this->heatTransferModels_,
353 heatTransferModelIter
354 )
355 {
356 const phasePair& pair
357 (
358 this->phasePairs_[heatTransferModelIter.key()]
359 );
360
361 const phaseModel& phase1 = pair.phase1();
362 const phaseModel& phase2 = pair.phase2();
363
364 const volScalarField& T1(phase1.thermo().T());
365 const volScalarField& T2(phase2.thermo().T());
366
367 const volScalarField& he1(phase1.thermo().he());
368 const volScalarField& he2(phase2.thermo().he());
369
370 const volScalarField& p(phase1.thermo().p());
371
372 volScalarField& iDmdt(*this->iDmdt_[pair]);
373 volScalarField& Tf(*this->Tf_[pair]);
374
375 const volScalarField Tsat(saturationModel_->Tsat(phase1.thermo().p()));
376
378 (
379 he1.member() == "e"
380 ? phase1.thermo().he(p, Tsat) + p/phase1.rho()
381 : phase1.thermo().he(p, Tsat)
382 );
384 (
385 he2.member() == "e"
386 ? phase2.thermo().he(p, Tsat) + p/phase2.rho()
387 : phase2.thermo().he(p, Tsat)
388 );
389
391 (
392 he1.member() == "e"
393 ? he1 + p/phase1.rho()
395 );
396
398 (
399 he2.member() == "e"
400 ? he2 + p/phase2.rho()
402 );
403
405 (
406 (neg0(iDmdt)*hf2 + pos(iDmdt)*h2)
407 - (pos0(iDmdt)*hf1 + neg(iDmdt)*h1)
408 );
409
410 volScalarField iDmdtNew(iDmdt);
411
412 if (phaseChange_)
413 {
414 volScalarField H1(heatTransferModelIter().first()->K(0));
415 volScalarField H2(heatTransferModelIter().second()->K(0));
416
417 iDmdtNew = (H1*(Tsat - T1) + H2*(Tsat - T2))/L;
418 }
419 else
420 {
421 iDmdtNew == dimensionedScalar(iDmdt.dimensions());
422 }
423
424 volScalarField H1(heatTransferModelIter().first()->K());
425 volScalarField H2(heatTransferModelIter().second()->K());
426
427 // Limit the H[12] to avoid /0
428 H1.max(SMALL);
429 H2.max(SMALL);
430
431 Tf = (H1*T1 + H2*T2 + iDmdtNew*L)/(H1 + H2);
432
433 Info<< "Tf." << pair.name()
434 << ": min = " << gMin(Tf.primitiveField())
435 << ", mean = " << gAverage(Tf.primitiveField())
436 << ", max = " << gMax(Tf.primitiveField())
437 << endl;
438
439 scalar iDmdtRelax(this->mesh().fieldRelaxationFactor("iDmdt"));
440 iDmdt = (1 - iDmdtRelax)*iDmdt + iDmdtRelax*iDmdtNew;
441
442 if (phaseChange_)
443 {
444 Info<< "iDmdt." << pair.name()
445 << ": min = " << gMin(iDmdt.primitiveField())
446 << ", mean = " << gAverage(iDmdt.primitiveField())
447 << ", max = " << gMax(iDmdt.primitiveField())
448 << ", integral = " << fvc::domainIntegrate(iDmdt).value()
449 << endl;
450 }
451
452 volScalarField& wDmdt(*this->wDmdt_[pair]);
453 volScalarField& wMDotL(*this->wMDotL_[pair]);
454 wDmdt *= 0.0;
455 wMDotL *= 0.0;
456
457 bool wallBoilingActive = false;
458
459 forAllConstIter(phasePair, pair, iter)
460 {
461 const phaseModel& phase = iter();
462 const phaseModel& otherPhase = iter.otherPhase();
463
464 if
465 (
466 phase.mesh().foundObject<volScalarField>
467 (
468 "alphat." + phase.name()
469 )
470 )
471 {
472 const volScalarField& alphat =
473 phase.mesh().lookupObject<volScalarField>
474 (
475 "alphat." + phase.name()
476 );
477
478 const fvPatchList& patches = this->mesh().boundary();
479 forAll(patches, patchi)
480 {
481 const fvPatch& currPatch = patches[patchi];
482
483 if
484 (
485 isA<alphatPhaseChangeWallFunction>
486 (
487 alphat.boundaryField()[patchi]
488 )
489 )
490 {
491 const alphatPhaseChangeWallFunction& PCpatch =
492 refCast<const alphatPhaseChangeWallFunction>
493 (
494 alphat.boundaryField()[patchi]
495 );
496
497 phasePairKey key(phase.name(), otherPhase.name());
498
499 if (PCpatch.activePhasePair(key))
500 {
501 wallBoilingActive = true;
502
503 const scalarField& patchDmdt =
504 PCpatch.dmdt(key);
505 const scalarField& patchMDotL =
506 PCpatch.mDotL(key);
507
508 const scalar sign
509 (
510 Pair<word>::compare(pair, key)
511 );
512
513 forAll(patchDmdt, facei)
514 {
515 const label faceCelli =
516 currPatch.faceCells()[facei];
517 wDmdt[faceCelli] -= sign*patchDmdt[facei];
518 wMDotL[faceCelli] -= sign*patchMDotL[facei];
519 }
520 }
521 }
522 }
523 }
524 }
525
526 if (wallBoilingActive)
527 {
528 Info<< "wDmdt." << pair.name()
529 << ": min = " << gMin(wDmdt.primitiveField())
530 << ", mean = " << gAverage(wDmdt.primitiveField())
531 << ", max = " << gMax(wDmdt.primitiveField())
532 << ", integral = " << fvc::domainIntegrate(wDmdt).value()
533 << endl;
534 }
535 }
536}
537
538
539template<class BasePhaseSystem>
541{
542 if (BasePhaseSystem::read())
543 {
544 bool readOK = true;
545
546 // Models ...
547
548 return readOK;
549 }
550 else
551 {
552 return false;
553 }
554}
555
556
557// ************************************************************************* //
CGAL::Exact_predicates_exact_constructions_kernel K
phaseModel & phase1
phaseModel & phase2
tmp< volScalarField > dmdt() const
Return the blended mass transfer rate.
const dimensionSet & dimensions() const
Return dimensions.
const Mesh & mesh() const
Return mesh.
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
void max(const dimensioned< Type > &dt)
Use the maximum of the field and specified value.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
bool set(const Key &key, T *ptr)
Assign a new entry, overwriting existing entries.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static word member(const word &name)
Return member (name without the extension)
Definition: IOobject.C:295
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:69
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
Class to provide interfacial heat and mass transfer between a number of phases according the interfac...
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
virtual void correctInterfaceThermo()
Correct the interface thermodynamics.
tmp< volScalarField > wDmdt(const phasePairKey &key) const
Return the boundary mass transfer rate for a pair.
wMDotLTable wMDotL_
Boundary thermal energy transfer rate.
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
iDmdtTable iDmdt_
Interfacial Mass transfer rate.
tmp< volScalarField > iDmdt(const phasePairKey &key) const
Return the interfacial mass transfer rate for a pair.
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
wDmdtTable wDmdt_
Boundary Mass transfer rate.
virtual bool read()
Read base phaseProperties dictionary.
const saturationModel & saturation() const
Return the saturationModel.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
virtual const volScalarField & T() const
Temperature [K].
Definition: basicThermo.C:622
virtual volScalarField & p()
Pressure [Pa].
Definition: basicThermo.C:610
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
Abstract base-class for all alphatWallFunctions supporting phase-change.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:712
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:113
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:61
const dimensionedScalar & rho() const
Definition: phaseModel.H:171
const phaseModel & otherPhase() const
Return the other phase in this two-phase system.
Definition: phaseModel.C:218
const word & name() const
Definition: phaseModel.H:144
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
An ordered or unorder pair of phase names. Typically specified as follows.
Definition: phasePairKey.H:68
bool ordered() const noexcept
Return the ordered flag.
Definition: phasePairKey.H:98
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:56
virtual word name() const
Pair name.
Definition: phasePair.C:69
const multiphaseInter::phaseModel & phase1() const
Definition: phasePairI.H:30
const multiphaseInter::phaseModel & phase2() const
Definition: phasePairI.H:36
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:57
const word & name() const
Definition: phase.H:111
Lookup type of boundary radiation properties.
Definition: lookup.H:66
virtual tmp< volScalarField > rho() const
Density [kg/m^3].
Definition: rhoThermo.C:224
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
const polyBoundaryMesh & patches
dynamicFvMesh & mesh
Volume integrate volField creating a volField.
Calculate the finiteVolume matrix for implicit and explicit sources.
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sign(const dimensionedScalar &ds)
const dimensionSet dimEnergy
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensionedScalar negPart(const dimensionedScalar &ds)
dimensionedScalar neg(const dimensionedScalar &ds)
Type gAverage(const FieldField< Field, Type > &f)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
dimensionedScalar neg0(const dimensionedScalar &ds)
Type gMin(const FieldField< Field, Type > &f)
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
const dimensionSet dimDensity
Type gMax(const FieldField< Field, Type > &f)
dimensionedScalar posPart(const dimensionedScalar &ds)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:381
volScalarField & he2
Definition: EEqns.H:3
const vector L(dict.get< vector >("L"))