HeatTransferPhaseSystem.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 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 
31 #include "HashPtrTable.H"
32 
33 #include "fvcDiv.H"
34 #include "fvmSup.H"
35 #include "fvMatrix.H"
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 template<class BasePhaseSystem>
41 (
42  const fvMesh& mesh
43 )
44 :
45  BasePhaseSystem(mesh)
46 {
47  this->generatePairsAndSubModels("heatTransfer", heatTransferModels_);
48 }
49 
50 
51 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
52 
53 template<class BasePhaseSystem>
55 {}
56 
57 
58 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
59 
60 template<class BasePhaseSystem>
62 (
63  const phaseModel& phase
64 ) const
65 {
66  return false;
67 }
68 
69 
70 template<class BasePhaseSystem>
73 (
74  const phasePairKey& key
75 ) const
76 {
77  return tmp<volScalarField>
78  (
79  new volScalarField
80  (
81  IOobject
82  (
84  (
85  "dmdt",
86  this->phasePairs_[key]->name()
87  ),
88  this->mesh().time().timeName(),
89  this->mesh().time()
90  ),
91  this->mesh(),
93  )
94  );
95 }
96 
97 
98 template<class BasePhaseSystem>
101 (
102  const Foam::phaseModel& phase
103 ) const
104 {
105  return tmp<volScalarField>(NULL);
106 }
107 
108 
109 template<class BasePhaseSystem>
112 {
113  autoPtr<phaseSystem::heatTransferTable> eqnsPtr
114  (
116  );
117 
118  phaseSystem::heatTransferTable& eqns = eqnsPtr();
119 
120  forAll(this->phaseModels_, phasei)
121  {
122  const phaseModel& phase = this->phaseModels_[phasei];
123 
124  eqns.insert
125  (
126  phase.name(),
127  new fvScalarMatrix(phase.thermo().he(), dimEnergy/dimTime)
128  );
129  }
130 
132  (
135  heatTransferModelIter
136  )
137  {
138  const volScalarField K(heatTransferModelIter()->K());
139 
140  const phasePair& pair(this->phasePairs_[heatTransferModelIter.key()]);
141 
142  const phaseModel* phase = &pair.phase1();
143  const phaseModel* otherPhase = &pair.phase2();
144 
145  forAllConstIter(phasePair, pair, iter)
146  {
147  const volScalarField& he(phase->thermo().he());
148  volScalarField Cpv(phase->thermo().Cpv());
149 
150  *eqns[phase->name()] +=
151  K*(otherPhase->thermo().T() - phase->thermo().T() + he/Cpv)
152  - fvm::Sp(K/Cpv, he);
153 
154  Swap(phase, otherPhase);
155  }
156  }
157 
158  return eqnsPtr;
159 }
160 
161 
162 template<class BasePhaseSystem>
165 {
166  autoPtr<phaseSystem::massTransferTable> eqnsPtr
167  (
169  );
170 
171  return eqnsPtr;
172 }
173 
174 
175 template<class BasePhaseSystem>
177 {
178  if (BasePhaseSystem::read())
179  {
180  bool readOK = true;
181 
182  // Models ...
183 
184  return readOK;
185  }
186  else
187  {
188  return false;
189  }
190 }
191 
192 
193 // ************************************************************************* //
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
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual autoPtr< phaseSystem::massTransferTable > massTransfer() const
Return the mass transfer matrices.
CGAL::Exact_predicates_exact_constructions_kernel K
virtual tmp< volScalarField > dmdt(const phasePairKey &key) const
Return the interfacial mass flow rate.
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
dynamicFvMesh & mesh
void Swap(T &a, T &b)
Definition: Swap.H:43
static word groupName(Name name, const word &group)
HashTable< HashTable< autoPtr< BlendedInterfacialModel< heatTransferModel > > >, phasePairKey, phasePairKey::hash > heatTransferModelTable
virtual ~HeatTransferPhaseSystem()
Destructor.
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Calculate the divergence of the given field.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
const dimensionSet dimDensity
virtual bool transfersMass(const phaseModel &phase) const
Return true if there is mass transfer for phase.
heatTransferModelTable heatTransferModels_
Heat transfer models.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
HashPtrTable< fvScalarMatrix, word, string::hash > massTransferTable
Definition: phaseSystem.H:118
A class for managing temporary objects.
Definition: PtrList.H:54
HeatTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
virtual bool read()
Read base phaseProperties dictionary.
Calculate the matrix for implicit and explicit sources.