PopulationBalancePhaseSystem.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) 2017-2018 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
29 
30 
31 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
32 
33 template<class BasePhaseSystem>
36 (
37  const phasePairKey& key
38 ) const
39 {
40  if (!pDmdt_.found(key))
41  {
42  return phaseSystem::dmdt(key);
43  }
44 
45  const scalar pDmdtSign(Pair<word>::compare(pDmdt_.find(key).key(), key));
46 
47  return pDmdtSign**pDmdt_[key];
48 }
49 
50 
51 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52 
53 template<class BasePhaseSystem>
56 (
57  const fvMesh& mesh
58 )
59 :
60  BasePhaseSystem(mesh),
61 
62  populationBalances_
63  (
64  this->lookup("populationBalances"),
66  )
67 {
68  forAll(populationBalances_, i)
69  {
71  populationBalances_[i];
72 
74  {
75  const phasePairKey& key = iter.key();
76 
77  if (!this->phasePairs_.found(key))
78  {
79  this->phasePairs_.insert
80  (
81  key,
83  (
84  new phasePair
85  (
86  this->phaseModels_[key.first()],
87  this->phaseModels_[key.second()]
88  )
89  )
90  );
91  }
92  }
93  }
94 
96  (
98  this->phasePairs_,
99  phasePairIter
100  )
101  {
102  const phasePair& pair(phasePairIter());
103 
104  if (pair.ordered())
105  {
106  continue;
107  }
108 
109  // Initially assume no mass transfer
110  pDmdt_.set
111  (
112  pair,
113  new volScalarField
114  (
115  IOobject
116  (
117  IOobject::groupName("pDmdt", pair.name()),
118  this->mesh().time().timeName(),
119  this->mesh(),
120  IOobject::READ_IF_PRESENT,
121  IOobject::AUTO_WRITE
122  ),
123  this->mesh(),
125  )
126  );
127  }
128 }
129 
130 
131 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
132 
133 template<class BasePhaseSystem>
136 {}
137 
138 
139 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
140 
141 template<class BasePhaseSystem>
144 (
145  const phasePairKey& key
146 ) const
147 {
148  return BasePhaseSystem::dmdt(key) + this->pDmdt(key);
149 }
150 
151 
152 template<class BasePhaseSystem>
155 {
156  PtrList<volScalarField> dmdts(BasePhaseSystem::dmdts());
157 
158  forAllConstIter(pDmdtTable, pDmdt_, pDmdtIter)
159  {
160  const phasePair& pair = this->phasePairs_[pDmdtIter.key()];
161  const volScalarField& pDmdt = *pDmdtIter();
162 
163  this->addField(pair.phase1(), "dmdt", pDmdt, dmdts);
164  this->addField(pair.phase2(), "dmdt", - pDmdt, dmdts);
165  }
166 
167  return dmdts;
168 }
169 
170 
171 template<class BasePhaseSystem>
174 {
177 
178  phaseSystem::massTransferTable& eqns = eqnsPtr();
179 
181  (
183  this->phasePairs_,
184  phasePairIter
185  )
186  {
187  const phasePair& pair(phasePairIter());
188 
189  if (pair.ordered())
190  {
191  continue;
192  }
193 
194  const phaseModel& phase = pair.phase1();
195  const phaseModel& otherPhase = pair.phase2();
196 
197  // Note that the phase YiEqn does not contain a continuity error term,
198  // so these additions represent the entire mass transfer
199 
200  const volScalarField dmdt(this->pDmdt(pair));
201  const volScalarField dmdt12(negPart(dmdt));
202  const volScalarField dmdt21(posPart(dmdt));
203 
204  const PtrList<volScalarField>& Yi = phase.Y();
205 
206  forAll(Yi, i)
207  {
208  const word name
209  (
210  IOobject::groupName(Yi[i].member(), phase.name())
211  );
212 
213  const word otherName
214  (
215  IOobject::groupName(Yi[i].member(), otherPhase.name())
216  );
217 
218  *eqns[name] +=
219  dmdt21*eqns[otherName]->psi()
220  + fvm::Sp(dmdt12, eqns[name]->psi());
221 
222  *eqns[otherName] -=
223  dmdt12*eqns[name]->psi()
224  + fvm::Sp(dmdt21, eqns[otherName]->psi());
225  }
226  }
227 
228  return eqnsPtr;
229 }
230 
231 
232 template<class BasePhaseSystem>
234 {
235  if (BasePhaseSystem::read())
236  {
237  bool readOK = true;
238 
239  // Models ...
240 
241  return readOK;
242  }
243  else
244  {
245  return false;
246  }
247 }
248 
249 
250 template<class BasePhaseSystem>
252 {
254 
255  forAll(populationBalances_, i)
256  {
257  populationBalances_[i].solve();
258  }
259 }
260 
261 
262 // ************************************************************************* //
PopulationBalancePhaseSystem.H
Foam::phaseModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:54
Foam::PopulationBalancePhaseSystem::massTransfer
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
Definition: PopulationBalancePhaseSystem.C:173
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::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::phase
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phase.H:54
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::dimDensity
const dimensionSet dimDensity
Foam::posPart
dimensionedScalar posPart(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:221
Sp
zeroField Sp
Definition: alphaSuSp.H:2
Foam::glTF::key
auto key(const Type &t) -> typename std::enable_if< std::is_enum< Type >::value, typename std::underlying_type< Type >::type >::type
Definition: foamGltfBase.H:108
Foam::diameterModels::populationBalanceModel::iNew
Return a pointer to a new populationBalanceModel object created on.
Definition: populationBalanceModel.H:357
solve
CEqn solve()
forAllConstIter
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:344
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::blockMeshTools::read
void read(Istream &, label &val, const dictionary &)
In-place read with dictionary lookup.
Definition: blockMeshTools.C:57
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::phasePairKey::ordered
bool ordered() const noexcept
Return the ordered flag.
Definition: phasePairKey.H:98
massTransfer
phaseSystem::massTransferTable & massTransfer(massTransferPtr())
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
Foam::phasePairKey
An ordered or unorder pair of phase names. Typically specified as follows.
Definition: phasePairKey.H:65
Foam::radiation::lookup
Lookup type of boundary radiation properties.
Definition: lookup.H:63
Foam::PopulationBalancePhaseSystem::dmdts
virtual PtrList< volScalarField > dmdts() const
Return the mass transfer rates for each phase.
Definition: PopulationBalancePhaseSystem.C:154
Foam::PopulationBalancePhaseSystem::dmdt
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the mass transfer rate for a pair.
Definition: PopulationBalancePhaseSystem.C:144
Foam::diameterModels::populationBalanceModel
Class that solves the univariate population balance equation by means of a class method (also called ...
Definition: populationBalanceModel.H:179
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::PopulationBalancePhaseSystem::solve
virtual void solve()
Solve all population balance equations.
Definition: PopulationBalancePhaseSystem.C:251
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam::phase::name
const word & name() const
Definition: phase.H:111
Foam::HashTable< autoPtr< phasePair >, phasePairKey, phasePairKey::hash >
Foam::phaseModel::name
const word & name() const
Definition: phaseModel.H:140
Foam::autoPtr< phasePair >
Foam::negPart
dimensionedScalar negPart(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:232
Foam::Pair
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: Pair.H:54
Foam::HashPtrTable< volScalarField, phasePairKey, phasePairKey::hash >
Foam::PopulationBalancePhaseSystem::pDmdt
virtual tmp< volScalarField > pDmdt(const phasePairKey &key) const
Return the population balance mass transfer rate.
Definition: PopulationBalancePhaseSystem.C:36
Foam::phasePair::name
virtual word name() const
Pair name.
Definition: phasePair.C:69
Foam::diameterModels::populationBalanceModel::phasePairs
const phasePairTable & phasePairs() const
Return list of unordered phasePairs in this populationBalance.
Definition: populationBalanceModelI.H:91
Foam::phasePair::phase2
const phaseModel & phase2() const
Definition: phasePairI.H:36
Foam::PopulationBalancePhaseSystem::read
virtual bool read()
Read base phaseProperties dictionary.
Definition: PopulationBalancePhaseSystem.C:233
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::PopulationBalancePhaseSystem::PopulationBalancePhaseSystem
PopulationBalancePhaseSystem(const fvMesh &)
Construct from fvMesh.
Definition: PopulationBalancePhaseSystem.C:56
Foam::phasePair::phase1
const phaseModel & phase1() const
Definition: phasePairI.H:30
Foam::GeometricField< scalar, fvPatchField, volMesh >
psi
const volScalarField & psi
Definition: createFieldRefs.H:1
Foam::PopulationBalancePhaseSystem::~PopulationBalancePhaseSystem
virtual ~PopulationBalancePhaseSystem()
Destructor.
Definition: PopulationBalancePhaseSystem.C:135