effectivenessHeatExchangerSource.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) 2013-2018 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 "fvMesh.H"
28 #include "fvMatrix.H"
30 #include "basicThermo.H"
31 #include "coupledPolyPatch.H"
32 #include "surfaceInterpolate.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace fv
39 {
40  defineTypeNameAndDebug(effectivenessHeatExchangerSource, 0);
42  (
43  option,
44  effectivenessHeatExchangerSource,
45  dictionary
46  );
47 }
48 }
49 
50 
51 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52 
53 void Foam::fv::effectivenessHeatExchangerSource::initialise()
54 {
55  const faceZone& fZone = mesh_.faceZones()[zoneID_];
56 
57  faceId_.setSize(fZone.size());
58  facePatchId_.setSize(fZone.size());
59  faceSign_.setSize(fZone.size());
60 
61  label count = 0;
62  forAll(fZone, i)
63  {
64  label facei = fZone[i];
65  label faceId = -1;
66  label facePatchId = -1;
67  if (mesh_.isInternalFace(facei))
68  {
69  faceId = facei;
70  facePatchId = -1;
71  }
72  else
73  {
74  facePatchId = mesh_.boundaryMesh().whichPatch(facei);
75  const polyPatch& pp = mesh_.boundaryMesh()[facePatchId];
76  if (isA<coupledPolyPatch>(pp))
77  {
78  if (refCast<const coupledPolyPatch>(pp).owner())
79  {
80  faceId = pp.whichFace(facei);
81  }
82  else
83  {
84  faceId = -1;
85  }
86  }
87  else if (!isA<emptyPolyPatch>(pp))
88  {
89  faceId = facei - pp.start();
90  }
91  else
92  {
93  faceId = -1;
94  facePatchId = -1;
95  }
96  }
97 
98  if (faceId >= 0)
99  {
100  if (fZone.flipMap()[i])
101  {
102  faceSign_[count] = -1;
103  }
104  else
105  {
106  faceSign_[count] = 1;
107  }
108  faceId_[count] = faceId;
109  facePatchId_[count] = facePatchId;
110  count++;
111  }
112  }
113  faceId_.setSize(count);
114  facePatchId_.setSize(count);
115  faceSign_.setSize(count);
116 
117  calculateTotalArea(faceZoneArea_);
118 }
119 
120 
121 void Foam::fv::effectivenessHeatExchangerSource::calculateTotalArea
122 (
123  scalar& area
124 )
125 {
126  area = 0;
127  forAll(faceId_, i)
128  {
129  label facei = faceId_[i];
130  if (facePatchId_[i] != -1)
131  {
132  label patchi = facePatchId_[i];
133  area += mesh_.magSf().boundaryField()[patchi][facei];
134  }
135  else
136  {
137  area += mesh_.magSf()[facei];
138  }
139  }
140  reduce(area, sumOp<scalar>());
141 }
142 
143 
144 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
145 
146 Foam::fv::effectivenessHeatExchangerSource::effectivenessHeatExchangerSource
147 (
148  const word& name,
149  const word& modelType,
150  const dictionary& dict,
151  const fvMesh& mesh
152 )
153 :
154  cellSetOption(name, modelType, dict, mesh),
155  secondaryMassFlowRate_(readScalar(coeffs_.lookup("secondaryMassFlowRate"))),
156  secondaryInletT_(readScalar(coeffs_.lookup("secondaryInletT"))),
157  primaryInletT_(readScalar(coeffs_.lookup("primaryInletT"))),
158  eTable_(),
159  UName_(coeffs_.lookupOrDefault<word>("U", "U")),
160  TName_(coeffs_.lookupOrDefault<word>("T", "T")),
161  phiName_(coeffs_.lookupOrDefault<word>("phi", "phi")),
162  faceZoneName_(coeffs_.lookup("faceZone")),
163  zoneID_(mesh_.faceZones().findZoneID(faceZoneName_)),
164  faceId_(),
165  facePatchId_(),
166  faceSign_(),
167  faceZoneArea_(0)
168 {
169  if (zoneID_ < 0)
170  {
172  << type() << " " << this->name() << ": "
173  << " Unknown face zone name: " << faceZoneName_
174  << ". Valid face zones are: " << mesh_.faceZones().names()
175  << nl << exit(FatalError);
176  }
177 
178  // Set the field name to that of the energy field from which the temperature
179  // is obtained
180 
181  const basicThermo& thermo =
182  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
183 
184  fieldNames_.setSize(1, thermo.he().name());
185 
186  applied_.setSize(1, false);
187 
188  eTable_.reset(new interpolation2DTable<scalar>(coeffs_));
189 
190  initialise();
191 }
192 
193 
194 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
195 
197 (
198  const volScalarField& rho,
199  fvMatrix<scalar>& eqn,
200  const label
201 )
202 {
203  const basicThermo& thermo =
204  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
205 
206  const surfaceScalarField Cpf(fvc::interpolate(thermo.Cp()));
207 
208  const surfaceScalarField& phi =
209  mesh_.lookupObject<surfaceScalarField>(phiName_);
210 
211  scalar totalphi = 0;
212  scalar CpfMean = 0;
213  forAll(faceId_, i)
214  {
215  label facei = faceId_[i];
216  if (facePatchId_[i] != -1)
217  {
218  label patchi = facePatchId_[i];
219  totalphi += phi.boundaryField()[patchi][facei]*faceSign_[i];
220 
221  CpfMean +=
222  Cpf.boundaryField()[patchi][facei]
223  *mesh_.magSf().boundaryField()[patchi][facei];
224  }
225  else
226  {
227  totalphi += phi[facei]*faceSign_[i];
228  CpfMean += Cpf[facei]*mesh_.magSf()[facei];
229  }
230  }
231  reduce(CpfMean, sumOp<scalar>());
232  reduce(totalphi, sumOp<scalar>());
233 
234  scalar Qt =
235  eTable_()(mag(totalphi), secondaryMassFlowRate_)
236  *(secondaryInletT_ - primaryInletT_)
237  *(CpfMean/faceZoneArea_)*mag(totalphi);
238 
239  const volScalarField& T = mesh_.lookupObject<volScalarField>(TName_);
240  const scalarField TCells(T, cells_);
241  scalar Tref = 0;
242  if (Qt > 0)
243  {
244  Tref = max(TCells);
245  reduce(Tref, maxOp<scalar>());
246  }
247  else
248  {
249  Tref = min(TCells);
250  reduce(Tref, minOp<scalar>());
251  }
252 
253  scalarField deltaTCells(cells_.size(), 0);
254  forAll(deltaTCells, i)
255  {
256  if (Qt > 0)
257  {
258  deltaTCells[i] = max(Tref - TCells[i], 0.0);
259  }
260  else
261  {
262  deltaTCells[i] = max(TCells[i] - Tref, 0.0);
263  }
264  }
265 
266  const volVectorField& U = mesh_.lookupObject<volVectorField>(UName_);
267  const scalarField& V = mesh_.V();
268  scalar sumWeight = 0;
269  forAll(cells_, i)
270  {
271  sumWeight += V[cells_[i]]*mag(U[cells_[i]])*deltaTCells[i];
272  }
273  reduce(sumWeight, sumOp<scalar>());
274 
275  if (this->V() > vSmall && mag(Qt) > vSmall)
276  {
277  scalarField& heSource = eqn.source();
278 
279  forAll(cells_, i)
280  {
281  heSource[cells_[i]] -=
282  Qt*V[cells_[i]]*mag(U[cells_[i]])*deltaTCells[i]/sumWeight;
283  }
284  }
285 
286  if (debug && Pstream::master())
287  {
288  Info<< indent << "Net mass flux [Kg/s] = " << totalphi << nl;
289  Info<< indent << "Total energy exchange [W] = " << Qt << nl;
290  Info<< indent << "Tref [K] = " << Tref << nl;
291  Info<< indent << "Efficiency : "
292  << eTable_()(mag(totalphi), secondaryMassFlowRate_) << endl;
293  }
294 }
295 
296 
298 {
299  if (cellSetOption::read(dict))
300  {
301  coeffs_.lookup("secondaryMassFlowRate") >> secondaryMassFlowRate_;
302  coeffs_.lookup("secondaryInletT") >> secondaryInletT_;
303  coeffs_.lookup("primaryInletT") >> primaryInletT_;
304 
305  return true;
306  }
307  else
308  {
309  return false;
310  }
311 }
312 
313 
314 // ************************************************************************* //
defineTypeNameAndDebug(option, 0)
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
label faceId(-1)
const word & name() const
Return name.
Definition: IOobject.H:297
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
surfaceScalarField & phi
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const Boundary & boundaryField() const
Return const-reference to the boundary field.
virtual bool read(const dictionary &dict)
Read dictionary.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:421
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
virtual bool read(const dictionary &dict)
Read source dictionary.
rhoReactionThermo & thermo
Definition: createFields.H:28
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
static const char nl
Definition: Ostream.H:265
Field< Type > & source()
Definition: fvMatrix.H:294
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
label patchi
U
Definition: pEqn.H:72
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
Cell-set options abtract base class. Provides a base set of controls, e.g.:
Definition: cellSetOption.H:69
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Scalar.
Namespace for OpenFOAM.