TwoResistanceHeatTransferPhaseSystem.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2015-2023 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 #include "heatTransferModel.H"
28 #include "fvmSup.H"
30 
31 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
32 
33 template<class BasePhaseSystem>
35 (
36  const phaseSystem::dmdtfTable& dmdtfs,
37  const phaseSystem::dmdtfTable& Tfs,
39  const latentHeatTransfer transfer,
41 ) const
42 {
44  (
45  dmdtfs,
46  Tfs,
47  scheme,
48  eqns
49  );
50 
51  // Loop the pairs
52  forAllConstIter(phaseSystem::dmdtfTable, dmdtfs, dmdtfIter)
53  {
54  const phaseInterface interface(*this, dmdtfIter.key());
55 
56  const volScalarField& dmdtf = *dmdtfIter();
57 
58  const volScalarField& Tf = *Tfs[interface];
59 
60  const phaseModel& phase1 = interface.phase1();
61  const phaseModel& phase2 = interface.phase2();
62  const rhoThermo& thermo1 = phase1.thermo();
63  const rhoThermo& thermo2 = phase2.thermo();
64 
65  // Transfer coefficients
67  heatTransferModels_[interface];
68  const volScalarField H1(heatTransferModel.modelInThe(phase1).K());
69  const volScalarField H2(heatTransferModel.modelInThe(phase2).K());
70  const volScalarField H1Fac(H1/(H1 + H2));
71  const volScalarField HEff(H1Fac*H2);
72 
73  // Latent heat contribution
74  switch (transfer)
75  {
76  case latentHeatTransfer::heat:
77  {
78  *eqns[phase1.name()] +=
79  - HEff*(thermo2.T() - thermo1.T()) + H1*(Tf - thermo1.T());
80 
81  *eqns[phase2.name()] +=
82  - HEff*(thermo1.T() - thermo2.T()) + H2*(Tf - thermo2.T());
83 
84  break;
85  }
86  case latentHeatTransfer::mass:
87  {
88  const volScalarField L(this->L(interface, dmdtf, Tf, scheme));
89 
90  *eqns[phase1.name()] += H1Fac*dmdtf*L;
91  *eqns[phase2.name()] += (1 - H1Fac)*dmdtf*L;
92 
93  break;
94  }
95  }
96  }
97 }
98 
99 
100 template<class BasePhaseSystem>
102 (
103  const phaseSystem::dmidtfTable& dmidtfs,
104  const phaseSystem::dmdtfTable& Tfs,
105  const latentHeatScheme scheme,
106  const latentHeatTransfer transfer,
108 ) const
109 {
111  (
112  dmidtfs,
113  Tfs,
114  scheme,
115  eqns
116  );
117 
118  // Loop the pairs
119  forAllConstIter(phaseSystem::dmidtfTable, dmidtfs, dmidtfIter)
120  {
121  const phaseInterface interface(*this, dmidtfIter.key());
122 
123  const volScalarField& Tf = *Tfs[interface];
124 
125  const phaseModel& phase1 = interface.phase1();
126  const phaseModel& phase2 = interface.phase2();
127  const rhoThermo& thermo1 = phase1.thermo();
128  const rhoThermo& thermo2 = phase2.thermo();
129 
130  // Transfer coefficients
132  heatTransferModels_[interface];
133  const volScalarField H1(heatTransferModel.modelInThe(phase1).K());
134  const volScalarField H2(heatTransferModel.modelInThe(phase2).K());
135  const volScalarField H1Fac(H1/(H1 + H2));
136  const volScalarField HEff(H1Fac*H2);
137 
138  // Loop the species
139  forAllConstIter(HashPtrTable<volScalarField>, *dmidtfIter(), dmidtfJter)
140  {
141  const word& specie = dmidtfJter.key();
142 
143  const volScalarField& dmidtf = *dmidtfJter();
144 
145  // Latent heat contribution
146  switch (transfer)
147  {
148  case latentHeatTransfer::heat:
149  {
150  // Do nothing. This term is handled outside the specie loop.
151 
152  break;
153  }
154  case latentHeatTransfer::mass:
155  {
156  const volScalarField Li
157  (
158  this->Li(interface, specie, dmidtf, Tf, scheme)
159  );
160 
161  *eqns[phase1.name()] += H1Fac*dmidtf*Li;
162  *eqns[phase2.name()] += (1 - H1Fac)*dmidtf*Li;
163 
164  break;
165  }
166  }
167  }
168 
169  // Latent heat contribution
170  switch (transfer)
171  {
172  case latentHeatTransfer::heat:
173  {
174  *eqns[phase1.name()] +=
175  - HEff*(thermo2.T() - thermo1.T()) + H1*(Tf - thermo1.T());
176 
177  *eqns[phase2.name()] +=
178  - HEff*(thermo1.T() - thermo2.T()) + H2*(Tf - thermo2.T());
179 
180  break;
181  }
182  case latentHeatTransfer::mass:
183  {
184  // Do nothing. This term is handled inside the specie loop.
185 
186  break;
187  }
188  }
189  }
190 }
191 
192 
193 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
194 
195 template<class BasePhaseSystem>
198 (
199  const fvMesh& mesh
200 )
201 :
202  HeatTransferPhaseSystem<BasePhaseSystem>(mesh)
203 {
204  this->generateInterfacialModels(heatTransferModels_);
205 
206  // Check that models have been specified on both sides of the interfaces
208  (
211  heatTransferModelIter
212  )
213  {
214  const phaseInterface& interface = heatTransferModelIter()->interface();
215 
216  forAllConstIter(phaseInterface, interface, iter)
217  {
218  if (!heatTransferModelIter()->haveModelInThe(iter()))
219  {
221  << "A heat transfer model for the " << iter().name()
222  << " side of the " << interface.name()
223  << " interface is not specified"
224  << exit(FatalError);
225  }
226  }
227  }
228 }
229 
230 
231 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
232 
233 template<class BasePhaseSystem>
236 {}
237 
238 
239 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
240 
241 template<class BasePhaseSystem>
244 heatTransfer() const
245 {
247  (
249  );
250 
251  phaseSystem::heatTransferTable& eqns = eqnsPtr();
252 
253  forAll(this->phaseModels_, phasei)
254  {
255  const phaseModel& phase = this->phaseModels_[phasei];
256 
257  eqns.insert
258  (
259  phase.name(),
260  new fvScalarMatrix(phase.thermo().he(), dimEnergy/dimTime)
261  );
262  }
263 
265  (
267  heatTransferModels_,
268  heatTransferModelIter
269  )
270  {
271  const sidedBlendedHeatTransferModel& model = heatTransferModelIter()();
272 
273  const phaseModel& phase1 = model.interface().phase1();
274  const phaseModel& phase2 = model.interface().phase2();
275  const volScalarField& he1 = phase1.thermo().he();
276  const volScalarField& he2 = phase2.thermo().he();
277  const volScalarField Cpv1(phase1.thermo().Cpv());
278  const volScalarField Cpv2(phase2.thermo().Cpv());
279 
280  const volScalarField H1(model.modelInThe(phase1).K());
281  const volScalarField H2(model.modelInThe(phase2).K());
282  const volScalarField HEff(H1*H2/(H1 + H2));
283 
284  *eqns[phase1.name()] +=
285  HEff*(phase2.thermo().T() - phase1.thermo().T())
286  + H1/Cpv1*he1 - fvm::Sp(H1/Cpv1, he1);
287  *eqns[phase2.name()] +=
288  HEff*(phase1.thermo().T() - phase2.thermo().T())
289  + H2/Cpv2*he2 - fvm::Sp(H2/Cpv2, he2);
290  }
291 
292  return eqnsPtr;
293 }
294 
295 
296 template<class BasePhaseSystem>
299 {
300  BasePhaseSystem::predictThermophysicalTransport();
301  correctInterfaceThermo();
302 }
303 
304 
305 template<class BasePhaseSystem>
308 {
309  BasePhaseSystem::correctThermophysicalTransport();
310 }
311 
312 
313 template<class BasePhaseSystem>
315 {
316  if (BasePhaseSystem::read())
317  {
318  bool readOK = true;
319 
320  // Models ...
321 
322  return readOK;
323  }
324  else
325  {
326  return false;
327  }
328 }
329 
330 
331 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Generic GeometricField class.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
const phaseInterface & interface() const
Access the interface.
const ModelType & modelInThe(const phaseModel &phase) const
Access the model within the given phase.
virtual autoPtr< phaseSystem::heatTransferTable > heatTransfer() const
Return the heat transfer matrices.
virtual void predictThermophysicalTransport()
Predict the energy transport e.g. alphat.
void addDmdtHefs(const phaseSystem::dmdtfTable &dmdtfs, const phaseSystem::dmdtfTable &Tfs, const latentHeatScheme scheme, const latentHeatTransfer transfer, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from bulk phase changes.
virtual void correctThermophysicalTransport()
Correct the energy transport e.g. alphat.
TwoResistanceHeatTransferPhaseSystem(const fvMesh &)
Construct from fvMesh.
heatTransferModelTable heatTransferModels_
Heat transfer models.
void addDmidtHefs(const phaseSystem::dmidtfTable &dmidtfs, const phaseSystem::dmdtfTable &Tfs, const latentHeatScheme scheme, const latentHeatTransfer transfer, phaseSystem::heatTransferTable &eqns) const
Add energy transfer terms which result from specie phase changes.
virtual bool read()
Read base phaseProperties dictionary.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
virtual const volScalarField & T() const =0
Temperature [K].
virtual const volScalarField & Cpv() const =0
Heat capacity at constant pressure/volume [J/kg/K].
virtual const volScalarField & he() const =0
Enthalpy/Internal energy [J/kg].
tmp< volScalarField > K() const
Return the heatTransfer coefficient K.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Model for heat transfer between phases.
tmp< volScalarField > K() const
The heat transfer function K used in the enthalpy equation.
latentHeatScheme
Enumeration for the latent heat scheme.
Class to represent an interface between phases. Derivations can further specify the configuration of ...
const phaseModel & phase1() const
Return phase 1.
const phaseModel & phase2() const
Return phase 2.
virtual const rhoThermo & thermo() const =0
Return the thermophysical model.
const word & name() const
Return the name of this phase.
Definition: phaseModel.C:145
Base-class for thermodynamic properties based on density.
Definition: rhoThermo.H:54
Base class of the thermophysical property types.
Definition: specie.H:66
latentHeatTransfer
Enumeration for the form of the latent heat transfer.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Calculate the matrix for implicit and explicit sources.
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
static tmp< surfaceInterpolationScheme< Type > > scheme(const surfaceScalarField &faceFlux, Istream &schemeData)
Return weighting factors for scheme given from Istream.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const VolField< Type > &)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const dimensionSet dimEnergy
const dimensionSet dimTime
error FatalError