externalWallHeatFluxTemperatureFvPatchScalarField.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) 2011-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 
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
31 
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  template<>
39  const char*
40  NamedEnum
41  <
43  3
44  >::names[] =
45  {
46  "power",
47  "flux",
48  "coefficient"
49  };
50 }
51 
52 const Foam::NamedEnum
53 <
55  3
57 
58 
59 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
60 
63 (
64  const fvPatch& p,
66 )
67 :
68  mixedFvPatchScalarField(p, iF),
69  temperatureCoupledBase(patch(), "undefined", "undefined", "undefined-K"),
70  mode_(fixedHeatFlux),
71  Q_(0),
72  Ta_(),
73  relaxation_(1),
74  emissivity_(0),
75  qrRelaxation_(1),
76  qrName_("undefined-qr"),
77  thicknessLayers_(),
78  kappaLayers_()
79 {
80  refValue() = 0;
81  refGrad() = 0;
82  valueFraction() = 1;
83 }
84 
85 
88 (
89  const fvPatch& p,
91  const dictionary& dict
92 )
93 :
94  mixedFvPatchScalarField(p, iF),
95  temperatureCoupledBase(patch(), dict),
96  mode_(operationModeNames.read(dict.lookup("mode"))),
97  Q_(0),
98  Ta_(),
99  relaxation_(dict.lookupOrDefault<scalar>("relaxation", 1)),
100  emissivity_(dict.lookupOrDefault<scalar>("emissivity", 0)),
101  qrRelaxation_(dict.lookupOrDefault<scalar>("qrRelaxation", 1)),
102  qrName_(dict.lookupOrDefault<word>("qr", "none")),
103  thicknessLayers_(),
104  kappaLayers_()
105 {
106  switch (mode_)
107  {
108  case fixedPower:
109  {
110  dict.lookup("Q") >> Q_;
111 
112  break;
113  }
114  case fixedHeatFlux:
115  {
116  q_ = scalarField("q", dict, p.size());
117 
118  break;
119  }
120  case fixedHeatTransferCoeff:
121  {
122  h_ = scalarField("h", dict, p.size());
123  Ta_ = Function1<scalar>::New("Ta", dict);
124 
125  if (dict.found("thicknessLayers"))
126  {
127  dict.lookup("thicknessLayers") >> thicknessLayers_;
128  dict.lookup("kappaLayers") >> kappaLayers_;
129  }
130 
131  break;
132  }
133  }
134 
135  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
136 
137  if (qrName_ != "none")
138  {
139  if (dict.found("qrPrevious"))
140  {
141  qrPrevious_ = scalarField("qrPrevious", dict, p.size());
142  }
143  else
144  {
145  qrPrevious_.setSize(p.size(), 0);
146  }
147  }
148 
149  if (dict.found("refValue"))
150  {
151  // Full restart
152  refValue() = scalarField("refValue", dict, p.size());
153  refGrad() = scalarField("refGradient", dict, p.size());
154  valueFraction() = scalarField("valueFraction", dict, p.size());
155  }
156  else
157  {
158  // Start from user entered data. Assume fixedValue.
159  refValue() = *this;
160  refGrad() = 0;
161  valueFraction() = 1;
162  }
163 }
164 
165 
168 (
170  const fvPatch& p,
172  const fvPatchFieldMapper& mapper
173 )
174 :
175  mixedFvPatchScalarField(ptf, p, iF, mapper),
176  temperatureCoupledBase(patch(), ptf),
177  mode_(ptf.mode_),
178  Q_(ptf.Q_),
179  Ta_(ptf.Ta_, false),
180  relaxation_(ptf.relaxation_),
181  emissivity_(ptf.emissivity_),
182  qrRelaxation_(ptf.qrRelaxation_),
183  qrName_(ptf.qrName_),
184  thicknessLayers_(ptf.thicknessLayers_),
185  kappaLayers_(ptf.kappaLayers_)
186 {
187  switch (mode_)
188  {
189  case fixedPower:
190  {
191  break;
192  }
193  case fixedHeatFlux:
194  {
195  q_.setSize(mapper.size());
196  q_.map(ptf.q_, mapper);
197 
198  break;
199  }
200  case fixedHeatTransferCoeff:
201  {
202  h_.setSize(mapper.size());
203  h_.map(ptf.h_, mapper);
204 
205  break;
206  }
207  }
208 
209  if (qrName_ != "none")
210  {
211  qrPrevious_.setSize(mapper.size());
212  qrPrevious_.map(ptf.qrPrevious_, mapper);
213  }
214 }
215 
216 
219 (
221 )
222 :
223  mixedFvPatchScalarField(tppsf),
224  temperatureCoupledBase(tppsf),
225  mode_(tppsf.mode_),
226  Q_(tppsf.Q_),
227  q_(tppsf.q_),
228  h_(tppsf.h_),
229  Ta_(tppsf.Ta_, false),
230  relaxation_(tppsf.relaxation_),
231  emissivity_(tppsf.emissivity_),
232  qrPrevious_(tppsf.qrPrevious_),
233  qrRelaxation_(tppsf.qrRelaxation_),
234  qrName_(tppsf.qrName_),
235  thicknessLayers_(tppsf.thicknessLayers_),
236  kappaLayers_(tppsf.kappaLayers_)
237 {}
238 
239 
242 (
245 )
246 :
247  mixedFvPatchScalarField(tppsf, iF),
248  temperatureCoupledBase(patch(), tppsf),
249  mode_(tppsf.mode_),
250  Q_(tppsf.Q_),
251  q_(tppsf.q_),
252  h_(tppsf.h_),
253  Ta_(tppsf.Ta_, false),
254  relaxation_(tppsf.relaxation_),
255  emissivity_(tppsf.emissivity_),
256  qrPrevious_(tppsf.qrPrevious_),
257  qrRelaxation_(tppsf.qrRelaxation_),
258  qrName_(tppsf.qrName_),
259  thicknessLayers_(tppsf.thicknessLayers_),
260  kappaLayers_(tppsf.kappaLayers_)
261 {}
262 
263 
264 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
265 
267 (
268  const fvPatchFieldMapper& m
269 )
270 {
271  mixedFvPatchScalarField::autoMap(m);
272 
273  switch (mode_)
274  {
275  case fixedPower:
276  {
277  break;
278  }
279  case fixedHeatFlux:
280  {
281  q_.autoMap(m);
282 
283  break;
284  }
285  case fixedHeatTransferCoeff:
286  {
287  h_.autoMap(m);
288 
289  break;
290  }
291  }
292 
293  if (qrName_ != "none")
294  {
295  qrPrevious_.autoMap(m);
296  }
297 }
298 
299 
301 (
302  const fvPatchScalarField& ptf,
303  const labelList& addr
304 )
305 {
306  mixedFvPatchScalarField::rmap(ptf, addr);
307 
309  refCast<const externalWallHeatFluxTemperatureFvPatchScalarField>(ptf);
310 
311  switch (mode_)
312  {
313  case fixedPower:
314  {
315  break;
316  }
317  case fixedHeatFlux:
318  {
319  q_.rmap(tiptf.q_, addr);
320 
321  break;
322  }
323  case fixedHeatTransferCoeff:
324  {
325  h_.rmap(tiptf.h_, addr);
326 
327  break;
328  }
329  }
330 
331  if (qrName_ != "none")
332  {
333  qrPrevious_.rmap(tiptf.qrPrevious_, addr);
334  }
335 }
336 
337 
339 {
340  if (updated())
341  {
342  return;
343  }
344 
345  const scalarField& Tp(*this);
346 
347  scalarField qr(Tp.size(), 0);
348  if (qrName_ != "none")
349  {
350  qr =
351  qrRelaxation_
352  *patch().lookupPatchField<volScalarField, scalar>(qrName_)
353  + (1 - qrRelaxation_)*qrPrevious_;
354 
355  qrPrevious_ = qr;
356  }
357 
358  switch (mode_)
359  {
360  case fixedPower:
361  {
362  refGrad() = (Q_/gSum(patch().magSf()) + qr)/kappa(Tp);
363  refValue() = 0;
364  valueFraction() = 0;
365 
366  break;
367  }
368  case fixedHeatFlux:
369  {
370  refGrad() = (q_ + qr)/kappa(Tp);
371  refValue() = 0;
372  valueFraction() = 0;
373 
374  break;
375  }
376  case fixedHeatTransferCoeff:
377  {
378  scalar totalSolidRes = 0;
379  if (thicknessLayers_.size())
380  {
381  forAll(thicknessLayers_, iLayer)
382  {
383  const scalar l = thicknessLayers_[iLayer];
384  if (kappaLayers_[iLayer] > 0)
385  {
386  totalSolidRes += l/kappaLayers_[iLayer];
387  }
388  }
389  }
390  scalarField hp(1/(1/h_ + totalSolidRes));
391 
392  const scalar Ta = Ta_->value(this->db().time().timeOutputValue());
393  scalarField hpTa(hp*Ta);
394 
395  if (emissivity_ > 0)
396  {
397  // Evaluate the radiative flux to the environment
398  // from the surface temperature ...
399  if (totalSolidRes > 0)
400  {
401  // ... including the effect of the solid wall thermal
402  // resistance
403  scalarField TpLambda(h_/(h_ + 1/totalSolidRes));
404  scalarField Ts(TpLambda*Tp + (1 - TpLambda)*Ta);
405  scalarField lambdaTa4(pow4((1 - TpLambda)*Ta));
406 
407  hp += emissivity_*sigma.value()*(pow4(Ts) - lambdaTa4)/Tp;
408  hpTa += emissivity_*sigma.value()*(lambdaTa4 + pow4(Ta));
409  }
410  else
411  {
412  // ... if there is no solid wall thermal resistance use
413  // the current wall temperature
414  hp += emissivity_*sigma.value()*pow3(Tp);
415  hpTa += emissivity_*sigma.value()*pow4(Ta);
416  }
417  }
418 
419  const scalarField kappaDeltaCoeffs
420  (
421  this->kappa(Tp)*patch().deltaCoeffs()
422  );
423 
424  refGrad() = 0;
425 
426  forAll(Tp, i)
427  {
428  if (qr[i] < 0)
429  {
430  const scalar hpmqr = hp[i] - qr[i]/Tp[i];
431 
432  refValue()[i] = hpTa[i]/hpmqr;
433  valueFraction()[i] = hpmqr/(hpmqr + kappaDeltaCoeffs[i]);
434  }
435  else
436  {
437  refValue()[i] = (hpTa[i] + qr[i])/hp[i];
438  valueFraction()[i] = hp[i]/(hp[i] + kappaDeltaCoeffs[i]);
439  }
440  }
441 
442  break;
443  }
444  }
445 
446  valueFraction() = relaxation_*valueFraction() + (1 - relaxation_);
447  refValue() = relaxation_*refValue() + (1 - relaxation_)*Tp;
448 
449  mixedFvPatchScalarField::updateCoeffs();
450 
451  if (debug)
452  {
453  const scalar Q = gSum(kappa(Tp)*patch().magSf()*snGrad());
454 
455  Info<< patch().boundaryMesh().mesh().name() << ':'
456  << patch().name() << ':'
457  << this->internalField().name() << " :"
458  << " heat transfer rate:" << Q
459  << " walltemperature "
460  << " min:" << gMin(*this)
461  << " max:" << gMax(*this)
462  << " avg:" << gAverage(*this)
463  << endl;
464  }
465 }
466 
467 
469 (
470  Ostream& os
471 ) const
472 {
474 
475  os.writeKeyword("mode")
476  << operationModeNames[mode_] << token::END_STATEMENT << nl;
478 
479  switch (mode_)
480  {
481  case fixedPower:
482  {
483  os.writeKeyword("Q")
484  << Q_ << token::END_STATEMENT << nl;
485 
486  break;
487  }
488  case fixedHeatFlux:
489  {
490  q_.writeEntry("q", os);
491 
492  break;
493  }
494  case fixedHeatTransferCoeff:
495  {
496  h_.writeEntry("h", os);
497  Ta_->writeData(os);
498 
499  if (relaxation_ < 1)
500  {
501  os.writeKeyword("relaxation")
502  << relaxation_ << token::END_STATEMENT << nl;
503  }
504 
505  if (emissivity_ > 0)
506  {
507  os.writeKeyword("emissivity")
508  << emissivity_ << token::END_STATEMENT << nl;
509  }
510 
511  if (thicknessLayers_.size())
512  {
513  thicknessLayers_.writeEntry("thicknessLayers", os);
514  kappaLayers_.writeEntry("kappaLayers", os);
515  }
516 
517  break;
518  }
519  }
520 
521  os.writeKeyword("qr")<< qrName_ << token::END_STATEMENT << nl;
522 
523  if (qrName_ != "none")
524  {
525  os.writeKeyword("qrRelaxation")
526  << qrRelaxation_ << token::END_STATEMENT << nl;
527 
528  qrPrevious_.writeEntry("qrPrevious", os);
529  }
530 
531  refValue().writeEntry("refValue", os);
532  refGrad().writeEntry("refGradient", os);
533  valueFraction().writeEntry("valueFraction", os);
534  writeEntry("value", os);
535 }
536 
537 
538 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
539 
540 namespace Foam
541 {
543  (
546  );
547 }
548 
549 // ************************************************************************* //
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:431
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Type gMin(const FieldField< Field, Type > &f)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual label size() const =0
Type gSum(const FieldField< Field, Type > &f)
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const Type & value() const
Return const reference to value.
virtual label size() const
Return size.
Definition: fvPatch.H:161
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
Type gMax(const FieldField< Field, Type > &f)
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
dimensionedScalar pow3(const dimensionedScalar &ds)
void setSize(const label)
Reset size of List.
Definition: List.C:281
Common functions used in temperature coupled boundaries.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Type gAverage(const FieldField< Field, Type > &f)
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:395
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
messageStream Info
externalWallHeatFluxTemperatureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
This boundary condition applies a heat flux condition to temperature on an external wall in one of th...
void write(Ostream &) const
Write.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
static autoPtr< Function1< Type > > New(const word &entryName, const dictionary &dict)
Selector.
Definition: Function1New.C:32
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576