thermalBaffle1DFvPatchScalarField.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-2019 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 
26 #include "volFields.H"
27 #include "surfaceFields.H"
29 #include "mapDistribute.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace compressible
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class solidType>
43 (
44  const fvPatch& p,
46 )
47 :
49  mixedFvPatchScalarField(p, iF),
50  TName_("T"),
51  baffleActivated_(true),
52  thickness_(p.size()),
53  Qs_(p.size()),
54  solidDict_(),
55  solidPtr_(nullptr),
56  qrPrevious_(p.size()),
57  qrRelaxation_(1),
58  qrName_("undefined-qr")
59 {}
60 
61 
62 template<class solidType>
65 (
67  const fvPatch& p,
69  const fvPatchFieldMapper& mapper
70 )
71 :
72  mappedPatchBase(p.patch(), ptf),
73  mixedFvPatchScalarField(ptf, p, iF, mapper),
74  TName_(ptf.TName_),
75  baffleActivated_(ptf.baffleActivated_),
76  thickness_(mapper(ptf.thickness_)),
77  Qs_(mapper(ptf.Qs_)),
78  solidDict_(ptf.solidDict_),
79  solidPtr_(ptf.solidPtr_),
80  qrPrevious_(mapper(ptf.qrPrevious_)),
81  qrRelaxation_(ptf.qrRelaxation_),
82  qrName_(ptf.qrName_)
83 {}
84 
85 
86 template<class solidType>
89 (
90  const fvPatch& p,
92  const dictionary& dict
93 )
94 :
95  mappedPatchBase(p.patch(), NEARESTPATCHFACE, dict),
96  mixedFvPatchScalarField(p, iF),
97  TName_("T"),
98  baffleActivated_(dict.lookupOrDefault<bool>("baffleActivated", true)),
99  thickness_(),
100  Qs_(p.size(), 0),
101  solidDict_(dict),
102  solidPtr_(),
103  qrPrevious_(p.size(), 0.0),
104  qrRelaxation_(dict.lookupOrDefault<scalar>("qrRelaxation", 1)),
105  qrName_(dict.lookupOrDefault<word>("qr", "none"))
106 {
107  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
108 
109  if (dict.found("thickness"))
110  {
111  thickness_ = scalarField("thickness", dict, p.size());
112  }
113 
114  if (dict.found("Qs"))
115  {
116  Qs_ = scalarField("Qs", dict, p.size());
117  }
118 
119  if (dict.found("qrPrevious"))
120  {
121  qrPrevious_ = scalarField("qrPrevious", dict, p.size());
122  }
123 
124  if (dict.found("refValue") && baffleActivated_)
125  {
126  // Full restart
127  refValue() = scalarField("refValue", dict, p.size());
128  refGrad() = scalarField("refGradient", dict, p.size());
129  valueFraction() = scalarField("valueFraction", dict, p.size());
130  }
131  else
132  {
133  // Start from user entered data. Assume zeroGradient.
134  refValue() = *this;
135  refGrad() = 0.0;
136  valueFraction() = 0.0;
137  }
138 
139 }
140 
141 
142 template<class solidType>
145 (
147 )
148 :
149  mappedPatchBase(ptf.patch().patch(), ptf),
150  mixedFvPatchScalarField(ptf),
151  TName_(ptf.TName_),
152  baffleActivated_(ptf.baffleActivated_),
153  thickness_(ptf.thickness_),
154  Qs_(ptf.Qs_),
155  solidDict_(ptf.solidDict_),
156  solidPtr_(ptf.solidPtr_),
157  qrPrevious_(ptf.qrPrevious_),
158  qrRelaxation_(ptf.qrRelaxation_),
159  qrName_(ptf.qrName_)
160 {}
161 
162 
163 template<class solidType>
166 (
169 )
170 :
171  mappedPatchBase(ptf.patch().patch(), ptf),
172  mixedFvPatchScalarField(ptf, iF),
173  TName_(ptf.TName_),
174  baffleActivated_(ptf.baffleActivated_),
175  thickness_(ptf.thickness_),
176  Qs_(ptf.Qs_),
177  solidDict_(ptf.solidDict_),
178  solidPtr_(ptf.solidPtr_),
179  qrPrevious_(ptf.qrPrevious_),
180  qrRelaxation_(ptf.qrRelaxation_),
181  qrName_(ptf.qrName_)
182 {}
183 
184 
185 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
186 
187 template<class solidType>
189 {
190  const label patchi = patch().index();
191 
192  const label nbrPatchi = samplePolyPatch().index();
193 
194  return (patchi < nbrPatchi);
195 }
196 
197 
198 template<class solidType>
200 {
201  if (this->owner())
202  {
203  if (solidPtr_.empty())
204  {
205  solidPtr_.reset(new solidType(solidDict_));
206  }
207  return solidPtr_();
208  }
209  else
210  {
211  const fvPatch& nbrPatch =
212  patch().boundaryMesh()[samplePolyPatch().index()];
213 
214  const thermalBaffle1DFvPatchScalarField& nbrField =
215  refCast<const thermalBaffle1DFvPatchScalarField>
216  (
217  nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
218  );
219 
220  return nbrField.solid();
221  }
222 }
223 
224 
225 template<class solidType>
227 baffleThickness() const
228 {
229  if (this->owner())
230  {
231  if (thickness_.size() != patch().size())
232  {
234  (
235  solidDict_
236  )<< " Field thickness has not been specified "
237  << " for patch " << this->patch().name()
238  << exit(FatalIOError);
239  }
240 
241  return thickness_;
242  }
243  else
244  {
245  const mapDistribute& mapDist = this->mappedPatchBase::map();
246 
247  const fvPatch& nbrPatch =
248  patch().boundaryMesh()[samplePolyPatch().index()];
249  const thermalBaffle1DFvPatchScalarField& nbrField =
250  refCast<const thermalBaffle1DFvPatchScalarField>
251  (
252  nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
253  );
254 
255  tmp<scalarField> tthickness
256  (
257  new scalarField(nbrField.baffleThickness())
258  );
259  scalarField& thickness = tthickness.ref();
260  mapDist.distribute(thickness);
261  return tthickness;
262  }
263 }
264 
265 
266 template<class solidType>
268 {
269  if (this->owner())
270  {
271  return Qs_;
272  }
273  else
274  {
275  const mapDistribute& mapDist = this->mappedPatchBase::map();
276 
277  const fvPatch& nbrPatch =
278  patch().boundaryMesh()[samplePolyPatch().index()];
279 
280  const thermalBaffle1DFvPatchScalarField& nbrField =
281  refCast<const thermalBaffle1DFvPatchScalarField>
282  (
283  nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
284  );
285 
286  tmp<scalarField> tQs(new scalarField(nbrField.Qs()));
287  scalarField& Qs = tQs.ref();
288  mapDist.distribute(Qs);
289  return tQs;
290  }
291 }
292 
293 
294 template<class solidType>
296 (
297  const fvPatchFieldMapper& m
298 )
299 {
301 
302  mixedFvPatchScalarField::autoMap(m);
303 
304  if (this->owner())
305  {
306  m(thickness_, thickness_);
307  m(Qs_, Qs_);
308  }
309 }
310 
311 
312 template<class solidType>
314 (
315  const fvPatchScalarField& ptf,
316  const labelList& addr
317 )
318 {
319  mixedFvPatchScalarField::rmap(ptf, addr);
320 
321  const thermalBaffle1DFvPatchScalarField& tiptf =
322  refCast<const thermalBaffle1DFvPatchScalarField>(ptf);
323 
324  if (this->owner())
325  {
326  thickness_.rmap(tiptf.thickness_, addr);
327  Qs_.rmap(tiptf.Qs_, addr);
328  }
329 }
330 
331 
332 template<class solidType>
334 {
335  if (updated())
336  {
337  return;
338  }
339  // Since we're inside initEvaluate/evaluate there might be processor
340  // comms underway. Change the tag we use.
341  int oldTag = UPstream::msgType();
342  UPstream::msgType() = oldTag+1;
343 
344  const mapDistribute& mapDist = this->mappedPatchBase::map();
345 
346  const label patchi = patch().index();
347 
348  const label nbrPatchi = samplePolyPatch().index();
349 
350  if (baffleActivated_)
351  {
352  const fvPatch& nbrPatch = patch().boundaryMesh()[nbrPatchi];
353 
354  const compressible::turbulenceModel& turbModel =
355  db().template lookupObject<compressible::turbulenceModel>
356  (
358  );
359 
360  // local properties
361  const scalarField kappaw(turbModel.kappaEff(patchi));
362 
363  const fvPatchScalarField& Tp =
364  patch().template lookupPatchField<volScalarField, scalar>(TName_);
365 
366 
367  scalarField qr(Tp.size(), 0.0);
368 
369  if (qrName_ != "none")
370  {
371  qr = patch().template lookupPatchField<volScalarField, scalar>
372  (qrName_);
373 
374  qr = qrRelaxation_*qr + (1.0 - qrRelaxation_)*qrPrevious_;
375  qrPrevious_ = qr;
376  }
377 
378  tmp<scalarField> Ti = patchInternalField();
379 
380  scalarField myKDelta(patch().deltaCoeffs()*kappaw);
381 
382  // nrb properties
383  scalarField nbrTp =
384  turbModel.transport().T().boundaryField()[nbrPatchi];
385  mapDist.distribute(nbrTp);
386 
387  // solid properties
388  scalarField kappas(patch().size(), 0.0);
389  forAll(kappas, i)
390  {
391  kappas[i] = solid().kappa(0.0, (Tp[i] + nbrTp[i])/2.0);
392  }
393 
394  scalarField KDeltaSolid(kappas/baffleThickness());
395 
396  scalarField alpha(KDeltaSolid - qr/Tp);
397 
398  valueFraction() = alpha/(alpha + myKDelta);
399 
400  refValue() = (KDeltaSolid*nbrTp + Qs()/2.0)/alpha;
401 
402  if (debug)
403  {
404  scalar Q = gAverage(kappaw*snGrad());
405  Info<< patch().boundaryMesh().mesh().name() << ':'
406  << patch().name() << ':'
407  << this->internalField().name() << " <- "
408  << nbrPatch.name() << ':'
409  << this->internalField().name() << " :"
410  << " heat[W]:" << Q
411  << " walltemperature "
412  << " min:" << gMin(*this)
413  << " max:" << gMax(*this)
414  << " avg:" << gAverage(*this)
415  << endl;
416  }
417  }
418 
419  // Restore tag
420  UPstream::msgType() = oldTag;
421 
422  mixedFvPatchScalarField::updateCoeffs();
423 }
424 
425 template<class solidType>
427 {
430 
431  if (this->owner())
432  {
433  writeEntry(os, "thickness", baffleThickness()());
434  writeEntry(os, "Qs", Qs()());
435  solid().write(os);
436  }
437 
438  writeEntry(os, "qrPrevious", qrPrevious_);
439  writeEntry(os, "qr", qrName_);
440  writeEntry(os, "qrRelaxation", qrRelaxation_);
441 }
442 
443 
444 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
445 
446 } // End namespace compressible
447 } // End namespace Foam
448 
449 // ************************************************************************* //
Foam::surfaceFields.
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:179
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
thermalBaffle1DFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
virtual tmp< volScalarField > kappaEff() const
Effective thermal turbulent diffusivity for temperature.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
Type gMin(const FieldField< Field, Type > &f)
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
virtual void write(Ostream &) const
Write as a dictionary.
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
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
const mapDistribute & map() const
Return reference to the parallel distribution map.
void write(Ostream &, const label, const dictionary &)
Write with dictionary lookup.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:137
Foam::fvPatchFieldMapper.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
virtual label size() const
Return size.
Definition: fvPatch.H:155
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Type gMax(const FieldField< Field, Type > &f)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
label patchi
Class containing processor-to-processor mapping information.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:331
Type gAverage(const FieldField< Field, Type > &f)
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:295
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
virtual const word & name() const
Return name.
Definition: fvPatch.H:143
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Namespace for OpenFOAM.
IOerror FatalIOError