phaseChange.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) 2018-2019 OpenFOAM Foundation
9  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 
29 #include "phaseChange.H"
31 #include "phaseSystem.H"
32 #include "phasePairKey.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace diameterModels
39 {
40 namespace driftModels
41 {
42  defineTypeNameAndDebug(phaseChange, 0);
43  addToRunTimeSelectionTable(driftModel, phaseChange, dictionary);
44 }
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
50 
52 (
53  const populationBalanceModel& popBal,
54  const dictionary& dict
55 )
56 :
57  driftModel(popBal, dict),
58  pairKeys_(dict.lookup("pairs")),
59  numberWeighted_(dict.getOrDefault<Switch>("numberWeighted", false)),
60  W_(pairKeys_.size())
61 {
62  const phaseSystem& fluid = popBal_.fluid();
63 
64  forAll(pairKeys_, i)
65  {
66  const phasePair& pair = fluid.phasePairs()[pairKeys_[i]];
67 
68  W_.set
69  (
70  i,
71  new volScalarField
72  (
73  IOobject
74  (
75  IOobject::groupName(type() + ":W", pair.name()),
76  popBal_.mesh().time().timeName(),
77  popBal_.mesh()
78  ),
79  popBal_.mesh(),
81  (
82  inv(numberWeighted_ ? dimVolume : dimLength),
83  Zero
84  )
85  )
86  );
87  }
88 }
89 
90 
91 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
92 
94 {
95  const phaseSystem& fluid = popBal_.fluid();
96 
97  forAll(pairKeys_, i)
98  {
99  W_[i] *= 0.0;
100  }
101 
102  forAll(pairKeys_, k)
103  {
104  if (fluid.phasePairs().found(pairKeys_[k]))
105  {
106  const phasePair& pair = fluid.phasePairs()[pairKeys_[k]];
107 
109  {
110  const velocityGroup& vgj = popBal_.velocityGroups()[j];
111  if (pair.contains(vgj.phase()))
112  {
113  forAll(vgj.sizeGroups(), i)
114  {
115  const sizeGroup& fi = vgj.sizeGroups()[i];
116 
117  W_[k] +=
118  fi*max(fi.phase(), SMALL)
119  /(numberWeighted_ ? fi.x() : fi.d());
120  }
121  }
122  }
123  }
124  }
125 }
126 
127 
129 (
130  volScalarField& driftRate,
131  const label i
132 )
133 {
134  const velocityGroup& vg = popBal_.sizeGroups()[i].VelocityGroup();
135 
136  forAll(pairKeys_, k)
137  {
138  const phasePair& pair =
139  popBal_.fluid().phasePairs()[pairKeys_[k]];
140 
141  if (pair.contains(vg.phase()))
142  {
143  const volScalarField& iDmdt =
144  popBal_.mesh().lookupObject<volScalarField>
145  (
146  IOobject::groupName("iDmdt", pair.name())
147  );
148 
149  const scalar iDmdtSign =
150  vg.phase().name() == pair.first() ? +1 : -1;
151 
152  const sizeGroup& fi = popBal_.sizeGroups()[i];
153 
154  tmp<volScalarField> dDriftRate
155  (
156  iDmdtSign*iDmdt/(fi.phase().rho()*W_[k])
157  );
158 
159  if (!numberWeighted_)
160  {
161  dDriftRate.ref() *= fi.x()/fi.d();
162  }
163 
164  driftRate += dDriftRate;
165  }
166  }
167 }
168 
169 
170 // ************************************************************************* //
Foam::diameterModels::populationBalanceModel::fluid
const phaseSystem & fluid() const
Return reference to the phaseSystem.
Definition: populationBalanceModelI.H:49
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::phasePair
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:53
Foam::Switch
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:77
Foam::diameterModels::sizeGroup::phase
const phaseModel & phase() const
Return const-reference to the phase.
Definition: sizeGroupI.H:38
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::diameterModels::driftModels::phaseChange::correct
virtual void correct()
Correct diameter independent expressions.
Definition: phaseChange.C:93
fluid
twoPhaseSystem & fluid
Definition: setRegionFluidFields.H:3
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::diameterModels::populationBalanceModel::velocityGroups
const UPtrList< velocityGroup > & velocityGroups() const
Return the velocityGroups belonging to this populationBalance.
Definition: populationBalanceModelI.H:77
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::inv
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
Definition: dimensionedSphericalTensor.C:73
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::diameterModels::velocityGroup
This diameterModel is intended for use with a populationBalanceModel in order to simulate polydispers...
Definition: velocityGroup.H:107
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::diameterModels::sizeGroup
This class represents a single sizeGroup belonging to a velocityGroup. The main property of a sizeGro...
Definition: sizeGroup.H:96
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::phasePair::contains
bool contains(const phaseModel &phase) const
Return true if this phasePair contains the given phase.
Definition: phasePairI.H:42
Foam::diameterModels::populationBalanceModel
Class that solves the univariate population balance equation by means of a class method (also called ...
Definition: populationBalanceModel.H:179
Foam::diameterModels::driftModels::defineTypeNameAndDebug
defineTypeNameAndDebug(constantDrift, 0)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::phaseModel::name
const word & name() const
Definition: phaseModel.H:140
Foam::diameterModels::sizeGroup::d
const dimensionedScalar & d() const
Return representative diameter of the sizeGroup.
Definition: sizeGroupI.H:52
phaseChange.H
Foam::phasePair::name
virtual word name() const
Pair name.
Definition: phasePair.C:69
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
Foam::roots::type
type
Types of root.
Definition: Roots.H:54
Foam::phaseSystem
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:66
Foam::phaseModel::rho
const dimensionedScalar & rho() const
Definition: phaseModel.H:167
Foam::diameterModels::driftModel
Base class for drift models.
Definition: driftModel.H:52
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::diameterModels::velocityGroup::sizeGroups
const PtrList< sizeGroup > & sizeGroups() const
Return sizeGroups belonging to this velocityGroup.
Definition: velocityGroupI.H:52
Foam::diameterModels::driftModel::popBal_
const populationBalanceModel & popBal_
Reference to the populationBalanceModel.
Definition: driftModel.H:59
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::diameterModels::sizeGroup::x
const dimensionedScalar & x() const
Return representative volume of the sizeGroup.
Definition: sizeGroupI.H:59
Foam::diameterModels::driftModels::phaseChange::phaseChange
phaseChange(const populationBalanceModel &popBal, const dictionary &dict)
Construct from a population balance model and a dictionary.
Definition: phaseChange.C:52
Foam::diameterModel::phase
const phaseModel & phase() const
Return the phase.
Definition: diameterModel.H:116
Foam::diameterModels::driftModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(driftModel, constantDrift, dictionary)
Foam::diameterModels::driftModels::phaseChange::addToDriftRate
virtual void addToDriftRate(volScalarField &driftRate, const label i)
Add to driftRate.
Definition: phaseChange.C:129