thermalBaffle1DFvPatchScalarField.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "volFields.H"
30#include "surfaceFields.H"
32#include "mapDistribute.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace compressible
39{
40
41// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42
43template<class solidType>
46(
47 const fvPatch& p,
49)
50:
51 mappedPatchBase(p.patch()),
52 mixedFvPatchScalarField(p, iF),
53 TName_("T"),
54 baffleActivated_(true),
55 thickness_(p.size()),
56 qs_(p.size()),
57 solidDict_(),
58 solidPtr_(nullptr),
59 qrPrevious_(p.size()),
60 qrRelaxation_(1),
61 qrName_("undefined-qr")
62{}
63
64
65template<class solidType>
68(
70 const fvPatch& p,
72 const fvPatchFieldMapper& mapper
73)
74:
75 mappedPatchBase(p.patch(), ptf),
76 mixedFvPatchScalarField(ptf, p, iF, mapper),
77 TName_(ptf.TName_),
78 baffleActivated_(ptf.baffleActivated_),
79 thickness_(ptf.thickness_, mapper),
80 qs_(ptf.qs_, mapper),
81 solidDict_(ptf.solidDict_),
82 solidPtr_(ptf.solidPtr_),
83 qrPrevious_(ptf.qrPrevious_, mapper),
84 qrRelaxation_(ptf.qrRelaxation_),
85 qrName_(ptf.qrName_)
86{}
87
88
89template<class solidType>
92(
93 const fvPatch& p,
95 const dictionary& dict
96)
97:
98 mappedPatchBase(p.patch(), NEARESTPATCHFACE, dict),
99 mixedFvPatchScalarField(p, iF),
100 TName_("T"),
101 baffleActivated_(dict.getOrDefault("baffleActivated", true)),
102 thickness_(),
103 qs_(p.size(), 0),
104 solidDict_(dict),
105 solidPtr_(),
106 qrPrevious_(p.size(), Zero),
107 qrRelaxation_
108 (
109 dict.getOrDefaultCompat("qrRelaxation", {{"relaxation", 1712}}, 1)
110 ),
111 qrName_(dict.getOrDefault<word>("qr", "none"))
112{
114
115 if (dict.found("thickness"))
116 {
117 thickness_ = scalarField("thickness", dict, p.size());
118 }
119
120 if (dict.found("qs"))
121 {
122 qs_ = scalarField("qs", dict, p.size());
123 }
124
125 if (dict.found("qrPrevious"))
126 {
127 qrPrevious_ = scalarField("qrPrevious", dict, p.size());
128 }
129
130 if (dict.found("refValue") && baffleActivated_)
131 {
132 // Full restart
133 refValue() = scalarField("refValue", dict, p.size());
134 refGrad() = scalarField("refGradient", dict, p.size());
135 valueFraction() = scalarField("valueFraction", dict, p.size());
136 }
137 else
138 {
139 // Start from user entered data. Assume zeroGradient.
140 refValue() = *this;
141 refGrad() = 0.0;
142 valueFraction() = 0.0;
143 }
144
145}
146
147
148template<class solidType>
151(
153)
154:
155 mappedPatchBase(ptf.patch().patch(), ptf),
156 mixedFvPatchScalarField(ptf),
157 TName_(ptf.TName_),
158 baffleActivated_(ptf.baffleActivated_),
159 thickness_(ptf.thickness_),
160 qs_(ptf.qs_),
161 solidDict_(ptf.solidDict_),
162 solidPtr_(ptf.solidPtr_),
163 qrPrevious_(ptf.qrPrevious_),
164 qrRelaxation_(ptf.qrRelaxation_),
165 qrName_(ptf.qrName_)
166{}
167
168
169template<class solidType>
172(
175)
176:
177 mappedPatchBase(ptf.patch().patch(), ptf),
178 mixedFvPatchScalarField(ptf, iF),
179 TName_(ptf.TName_),
180 baffleActivated_(ptf.baffleActivated_),
181 thickness_(ptf.thickness_),
182 qs_(ptf.qs_),
183 solidDict_(ptf.solidDict_),
184 solidPtr_(ptf.solidPtr_),
185 qrPrevious_(ptf.qrPrevious_),
186 qrRelaxation_(ptf.qrRelaxation_),
187 qrName_(ptf.qrName_)
188{}
189
190
191// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
192
193template<class solidType>
195{
196 const label patchi = patch().index();
197
198 const label nbrPatchi = samplePolyPatch().index();
199
200 return (patchi < nbrPatchi);
201}
202
203
204template<class solidType>
205const solidType& thermalBaffle1DFvPatchScalarField<solidType>::solid() const
206{
207 if (this->owner())
208 {
209 if (!solidPtr_)
210 {
211 solidPtr_.reset(new solidType(solidDict_));
212 }
213 return *solidPtr_;
214 }
215 else
216 {
217 const fvPatch& nbrPatch =
218 patch().boundaryMesh()[samplePolyPatch().index()];
219
220 const thermalBaffle1DFvPatchScalarField& nbrField =
221 refCast<const thermalBaffle1DFvPatchScalarField>
222 (
223 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
224 );
225
226 return nbrField.solid();
227 }
228}
229
230
231template<class solidType>
232tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::
233baffleThickness() const
234{
235 if (this->owner())
236 {
237 if (thickness_.size() != patch().size())
238 {
239 FatalIOErrorInFunction(solidDict_)
240 << "Field thickness has not been specified"
241 " for patch " << this->patch().name()
242 << exit(FatalIOError);
243 }
244
245 return thickness_;
246 }
247 else
248 {
249 const mapDistribute& mapDist = this->mappedPatchBase::map();
250
251 const fvPatch& nbrPatch =
252 patch().boundaryMesh()[samplePolyPatch().index()];
253 const thermalBaffle1DFvPatchScalarField& nbrField =
254 refCast<const thermalBaffle1DFvPatchScalarField>
255 (
256 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
257 );
258
259 tmp<scalarField> tthickness
260 (
261 new scalarField(nbrField.baffleThickness())
262 );
263 scalarField& thickness = tthickness.ref();
264 mapDist.distribute(thickness);
265 return tthickness;
266 }
267}
268
269
270template<class solidType>
271tmp<scalarField> thermalBaffle1DFvPatchScalarField<solidType>::qs() const
272{
273 if (this->owner())
274 {
275 return qs_;
276 }
277 else
278 {
279 const mapDistribute& mapDist = this->mappedPatchBase::map();
280
281 const fvPatch& nbrPatch =
282 patch().boundaryMesh()[samplePolyPatch().index()];
283
284 const thermalBaffle1DFvPatchScalarField& nbrField =
285 refCast<const thermalBaffle1DFvPatchScalarField>
286 (
287 nbrPatch.template lookupPatchField<volScalarField, scalar>(TName_)
288 );
289
290 tmp<scalarField> tqs(new scalarField(nbrField.qs()));
291 scalarField& qs = tqs.ref();
292 mapDist.distribute(qs);
293 return tqs;
294 }
295}
296
297
298template<class solidType>
300(
301 const fvPatchFieldMapper& m
302)
303{
305
306 mixedFvPatchScalarField::autoMap(m);
307
308 if (this->owner())
309 {
310 thickness_.autoMap(m);
311 qs_.autoMap(m);
312 }
313}
314
315
316template<class solidType>
318(
319 const fvPatchScalarField& ptf,
320 const labelList& addr
321)
322{
323 mixedFvPatchScalarField::rmap(ptf, addr);
324
326 refCast<const thermalBaffle1DFvPatchScalarField>(ptf);
327
328 if (this->owner())
329 {
330 thickness_.rmap(tiptf.thickness_, addr);
331 qs_.rmap(tiptf.qs_, addr);
332 }
333}
334
335
336template<class solidType>
338{
339 if (updated())
340 {
341 return;
342 }
343 // Since we're inside initEvaluate/evaluate there might be processor
344 // comms underway. Change the tag we use.
345 int oldTag = UPstream::msgType();
346 UPstream::msgType() = oldTag+1;
347
348 const mapDistribute& mapDist = this->mappedPatchBase::map();
349
350 const label patchi = patch().index();
351
352 const label nbrPatchi = samplePolyPatch().index();
353
354 if (baffleActivated_)
355 {
356 const fvPatch& nbrPatch = patch().boundaryMesh()[nbrPatchi];
357
358 const compressible::turbulenceModel& turbModel =
359 db().template lookupObject<compressible::turbulenceModel>
360 (
362 );
363
364 // local properties
365 const scalarField kappaw(turbModel.kappaEff(patchi));
366
367 const fvPatchScalarField& Tp =
368 patch().template lookupPatchField<volScalarField, scalar>(TName_);
369
370
371 scalarField qr(Tp.size(), Zero);
372
373 if (qrName_ != "none")
374 {
375 qr = patch().template lookupPatchField<volScalarField, scalar>
376 (qrName_);
377
378 qr = qrRelaxation_*qr + (1.0 - qrRelaxation_)*qrPrevious_;
379 qrPrevious_ = qr;
380 }
381
382 tmp<scalarField> Ti = patchInternalField();
383
384 scalarField myKDelta(patch().deltaCoeffs()*kappaw);
385
386 // nrb properties
387 scalarField nbrTp =
388 turbModel.transport().T().boundaryField()[nbrPatchi];
389 mapDist.distribute(nbrTp);
390
391 // solid properties
392 scalarField kappas(patch().size(), Zero);
393 forAll(kappas, i)
394 {
395 kappas[i] = solid().kappa(0.0, (Tp[i] + nbrTp[i])/2.0);
396 }
397
398 scalarField KDeltaSolid(kappas/baffleThickness());
399
400 scalarField alpha(KDeltaSolid - qr/Tp);
401
402 valueFraction() = alpha/(alpha + myKDelta);
403
404 refValue() = (KDeltaSolid*nbrTp + qs()/2.0)/alpha;
405
406 if (debug)
407 {
408 scalar Q = gAverage(kappaw*snGrad());
409 Info<< patch().boundaryMesh().mesh().name() << ':'
410 << patch().name() << ':'
411 << this->internalField().name() << " <- "
412 << nbrPatch.name() << ':'
413 << this->internalField().name() << " :"
414 << " heat[W]:" << Q
415 << " walltemperature "
416 << " min:" << gMin(*this)
417 << " max:" << gMax(*this)
418 << " avg:" << gAverage(*this)
419 << endl;
420 }
421 }
422
423 // Restore tag
424 UPstream::msgType() = oldTag;
425
426 mixedFvPatchScalarField::updateCoeffs();
427}
428
429template<class solidType>
431{
432 mixedFvPatchScalarField::write(os);
434
435 if (this->owner())
436 {
437 baffleThickness()().writeEntry("thickness", os);
438 qs()().writeEntry("qs", os);
439 solid().write(os);
440 }
441
442 qrPrevious_.writeEntry("qrPrevious", os);
443 os.writeEntry("qr", qrName_);
444 os.writeEntry("qrRelaxation", qrRelaxation_);
445}
446
447
448// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
449
450} // End namespace compressible
451} // End namespace Foam
452
453// ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual tmp< volScalarField > kappaEff() const
Return the effective turbulent thermal diffusivity for temperature.
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
static int & msgType() noexcept
Message tag of standard messages.
Definition: UPstream.H:556
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual void operator=(const UList< scalar > &)
Definition: fvPatchField.C:408
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
virtual const word & name() const
Return name.
Definition: fvPatch.H:173
Class containing processor-to-processor mapping information.
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
Determines a mapping between patch face centres and mesh cell or face centres and processors they're ...
const mapDistribute & map() const
Return reference to the parallel distribution map.
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:290
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
OBJstream os(runTime.globalPath()/outputName)
bool compressible
Definition: pEqn.H:2
const std::string patch
OpenFOAM patch number as a std::string.
Namespace for OpenFOAM.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
IOerror FatalIOError
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Type gMin(const FieldField< Field, Type > &f)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Type gMax(const FieldField< Field, Type > &f)
volScalarField & alpha
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.