effectivenessHeatExchangerSource.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) 2013-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 #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>("UName", "U")),
160  TName_(coeffs_.lookupOrDefault<word>("TName", "T")),
161  phiName_(coeffs_.lookupOrDefault<word>("phiName", "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  (
173  "effectivenessHeatExchangerSource::effectivenessHeatExchangerSource"
174  "("
175  "const word&, "
176  "const word&, "
177  "const dictionary&, "
178  "const fvMesh&"
179  ")"
180  )
181  << type() << " " << this->name() << ": "
182  << " Unknown face zone name: " << faceZoneName_
183  << ". Valid face zones are: " << mesh_.faceZones().names()
184  << nl << exit(FatalError);
185  }
186 
187  // Set the field name to that of the energy field from which the temperature
188  // is obtained
189 
190  const basicThermo& thermo =
191  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
192 
193  fieldNames_.setSize(1, thermo.he().name());
194 
195  applied_.setSize(1, false);
196 
197  eTable_.reset(new interpolation2DTable<scalar>(coeffs_));
198 
199  initialise();
200 }
201 
202 
203 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
204 
206 (
207  const volScalarField& rho,
208  fvMatrix<scalar>& eqn,
209  const label
210 )
211 {
212  const basicThermo& thermo =
213  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
214 
215  const surfaceScalarField Cpf(fvc::interpolate(thermo.Cp()));
216 
217  const surfaceScalarField& phi =
218  mesh_.lookupObject<surfaceScalarField>(phiName_);
219 
220  scalar totalphi = 0;
221  scalar CpfMean = 0;
222  forAll(faceId_, i)
223  {
224  label faceI = faceId_[i];
225  if (facePatchId_[i] != -1)
226  {
227  label patchI = facePatchId_[i];
228  totalphi += phi.boundaryField()[patchI][faceI]*faceSign_[i];
229 
230  CpfMean +=
231  Cpf.boundaryField()[patchI][faceI]
232  *mesh_.magSf().boundaryField()[patchI][faceI];
233  }
234  else
235  {
236  totalphi += phi[faceI]*faceSign_[i];
237  CpfMean += Cpf[faceI]*mesh_.magSf()[faceI];
238  }
239  }
240  reduce(CpfMean, sumOp<scalar>());
241  reduce(totalphi, sumOp<scalar>());
242 
243  scalar Qt =
244  eTable_()(mag(totalphi), secondaryMassFlowRate_)
245  *(secondaryInletT_ - primaryInletT_)
246  *(CpfMean/faceZoneArea_)*mag(totalphi);
247 
248  const volScalarField& T = mesh_.lookupObject<volScalarField>(TName_);
249  const scalarField TCells(T, cells_);
250  scalar Tref = 0;
251  if (Qt > 0)
252  {
253  Tref = max(TCells);
254  reduce(Tref, maxOp<scalar>());
255  }
256  else
257  {
258  Tref = min(TCells);
259  reduce(Tref, minOp<scalar>());
260  }
261 
262  scalarField deltaTCells(cells_.size(), 0);
263  forAll(deltaTCells, i)
264  {
265  if (Qt > 0)
266  {
267  deltaTCells[i] = max(Tref - TCells[i], 0.0);
268  }
269  else
270  {
271  deltaTCells[i] = max(TCells[i] - Tref, 0.0);
272  }
273  }
274 
275  const volVectorField& U = mesh_.lookupObject<volVectorField>(UName_);
276  const scalarField& V = mesh_.V();
277  scalar sumWeight = 0;
278  forAll(cells_, i)
279  {
280  sumWeight += V[cells_[i]]*mag(U[cells_[i]])*deltaTCells[i];
281  }
282  reduce(sumWeight, sumOp<scalar>());
283 
284  if (this->V() > VSMALL && mag(Qt) > VSMALL)
285  {
286  scalarField& heSource = eqn.source();
287 
288  forAll(cells_, i)
289  {
290  heSource[cells_[i]] -=
291  Qt*V[cells_[i]]*mag(U[cells_[i]])*deltaTCells[i]/sumWeight;
292  }
293  }
294 
295  if (debug && Pstream::master())
296  {
297  Info<< indent << "Net mass flux [Kg/s] = " << totalphi << nl;
298  Info<< indent << "Total energy exchange [W] = " << Qt << nl;
299  Info<< indent << "Tref [K] = " << Tref << nl;
300  Info<< indent << "Efficiency : "
301  << eTable_()(mag(totalphi), secondaryMassFlowRate_) << endl;
302  }
303 }
304 
305 
307 {
308  if (cellSetOption::read(dict))
309  {
310  coeffs_.lookup("secondaryMassFlowRate") >> secondaryMassFlowRate_;
311  coeffs_.lookup("secondaryInletT") >> secondaryInletT_;
312  coeffs_.lookup("primaryInletT") >> primaryInletT_;
313 
314  return true;
315  }
316  else
317  {
318  return false;
319  }
320 }
321 
322 
323 // ************************************************************************* //
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &vf, const surfaceScalarField &faceFlux, Istream &schemeData)
defineTypeNameAndDebug(cellSetOption, 0)
Cell-set options abtract base class. Provides a base set of controls, e.g.
Definition: cellSetOption.H:71
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual bool read(const dictionary &dict)
Read source dictionary.
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
A class for handling words, derived from string.
Definition: word.H:59
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Surface Interpolation.
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:68
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
Namespace for OpenFOAM.
virtual bool read(const dictionary &dict)
Read dictionary.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
label faceId(-1)
#define forAll(list, i)
Definition: UList.H:421
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
Macros for easy insertion into run-time selection tables.
const word & name() const
Return name.
Definition: IOobject.H:260
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldI)
Scalar.
labelList fv(nPoints)
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:398
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
Field< Type > & source()
Definition: fvMatrix.H:291
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:82
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].