turbulentTemperatureRadCoupledMixedFvPatchScalarField.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-2021 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"
30 #include "mappedPatchBase.H"
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace compressible
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
43 (
44  const fvPatch& p,
46 )
47 :
48  mixedFvPatchScalarField(p, iF),
49  temperatureCoupledBase(patch()),
50  TnbrName_("undefined-Tnbr"),
51  qrNbrName_("undefined-qrNbr"),
52  qrName_("undefined-qr"),
53  thicknessLayers_(0),
54  kappaLayers_(0),
55  contactRes_(0)
56 {
57  this->refValue() = 0.0;
58  this->refGrad() = 0.0;
59  this->valueFraction() = 1.0;
60 }
61 
62 
65 (
66  const fvPatch& p,
68  const dictionary& dict
69 )
70 :
71  mixedFvPatchScalarField(p, iF),
72  temperatureCoupledBase(patch(), dict),
73  TnbrName_(dict.lookupOrDefault<word>("Tnbr", "T")),
74  qrNbrName_(dict.lookupOrDefault<word>("qrNbr", "none")),
75  qrName_(dict.lookupOrDefault<word>("qr", "none")),
76  thicknessLayers_(0),
77  kappaLayers_(0),
78  contactRes_(0.0)
79 {
80  if (!isA<mappedPatchBase>(this->patch().patch()))
81  {
83  << "' not type '" << mappedPatchBase::typeName << "'"
84  << "\n for patch " << p.name()
85  << " of field " << internalField().name()
86  << " in file " << internalField().objectPath()
87  << exit(FatalError);
88  }
89 
90  if (dict.found("thicknessLayers"))
91  {
92  dict.lookup("thicknessLayers") >> thicknessLayers_;
93  dict.lookup("kappaLayers") >> kappaLayers_;
94 
95  if (thicknessLayers_.size() > 0)
96  {
97  // Calculate effective thermal resistance by harmonic averaging
98  forAll(thicknessLayers_, iLayer)
99  {
100  contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
101  }
102  contactRes_ = 1.0/contactRes_;
103  }
104  }
105 
106  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
107 
108  if (dict.found("refValue"))
109  {
110  // Full restart
111  refValue() = scalarField("refValue", dict, p.size());
112  refGrad() = scalarField("refGradient", dict, p.size());
113  valueFraction() = scalarField("valueFraction", dict, p.size());
114  }
115  else
116  {
117  // Start from user entered data. Assume fixedValue.
118  refValue() = *this;
119  refGrad() = 0.0;
120  valueFraction() = 1.0;
121  }
122 }
123 
124 
127 (
129  const fvPatch& p,
131  const fvPatchFieldMapper& mapper
132 )
133 :
134  mixedFvPatchScalarField(psf, p, iF, mapper),
135  temperatureCoupledBase(patch(), psf),
136  TnbrName_(psf.TnbrName_),
137  qrNbrName_(psf.qrNbrName_),
138  qrName_(psf.qrName_),
139  thicknessLayers_(psf.thicknessLayers_),
140  kappaLayers_(psf.kappaLayers_),
141  contactRes_(psf.contactRes_)
142 {}
143 
144 
147 (
150 )
151 :
152  mixedFvPatchScalarField(psf, iF),
153  temperatureCoupledBase(patch(), psf),
154  TnbrName_(psf.TnbrName_),
155  qrNbrName_(psf.qrNbrName_),
156  qrName_(psf.qrName_),
157  thicknessLayers_(psf.thicknessLayers_),
158  kappaLayers_(psf.kappaLayers_),
159  contactRes_(psf.contactRes_)
160 {}
161 
162 
163 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
164 
166 {
167  if (updated())
168  {
169  return;
170  }
171 
172  // Since we're inside initEvaluate/evaluate there might be processor
173  // comms underway. Change the tag we use.
174  int oldTag = UPstream::msgType();
175  UPstream::msgType() = oldTag+1;
176 
177  // Get the coupling information from the mappedPatchBase
178  const mappedPatchBase& mpp =
179  refCast<const mappedPatchBase>(patch().patch());
180  const polyMesh& nbrMesh = mpp.sampleMesh();
181  const label samplePatchi = mpp.samplePolyPatch().index();
182  const fvPatch& nbrPatch =
183  refCast<const fvMesh>(nbrMesh).boundary()[samplePatchi];
184 
185 
186  scalarField Tc(patchInternalField());
187  scalarField& Tp = *this;
188 
190 
191  const fvPatchScalarField& nbrTp =
192  nbrPatch.lookupPatchField<volScalarField, scalar>(TnbrName_);
193 
194  if (!isA<thisType>(nbrTp))
195  {
197  << "Patch field for " << internalField().name() << " on "
198  << patch().name() << " is of type " << thisType::typeName
199  << endl << "The neighbouring patch field " << TnbrName_ << " on "
200  << nbrPatch.name() << " is required to be the same, but is "
201  << "currently of type " << nbrTp.type() << exit(FatalError);
202  }
203 
204  const thisType& nbrField = refCast<const thisType>(nbrTp);
205 
206  // Swap to obtain full local values of neighbour internal field
207  scalarField TcNbr(nbrField.patchInternalField());
208  mpp.distribute(TcNbr);
209 
210 
211  // Swap to obtain full local values of neighbour K*delta
212  scalarField KDeltaNbr;
213  if (contactRes_ == 0.0)
214  {
215  KDeltaNbr = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
216  }
217  else
218  {
219  KDeltaNbr.setSize(nbrField.size(), contactRes_);
220  }
221  mpp.distribute(KDeltaNbr);
222 
223  scalarField KDelta(kappa(*this)*patch().deltaCoeffs());
224 
225  scalarField qr(Tp.size(), 0.0);
226  if (qrName_ != "none")
227  {
228  qr = patch().lookupPatchField<volScalarField, scalar>(qrName_);
229  }
230 
231  scalarField qrNbr(Tp.size(), 0.0);
232  if (qrNbrName_ != "none")
233  {
234  qrNbr = nbrPatch.lookupPatchField<volScalarField, scalar>(qrNbrName_);
235  mpp.distribute(qrNbr);
236  }
237 
238  valueFraction() = KDeltaNbr/(KDeltaNbr + KDelta);
239  refValue() = TcNbr;
240  refGrad() = (qr + qrNbr)/kappa(*this);
241 
242  mixedFvPatchScalarField::updateCoeffs();
243 
244  if (debug)
245  {
246  scalar Q = gSum(kappa(*this)*patch().magSf()*snGrad());
247 
248  Info<< patch().boundaryMesh().mesh().name() << ':'
249  << patch().name() << ':'
250  << this->internalField().name() << " <- "
251  << nbrMesh.name() << ':'
252  << nbrPatch.name() << ':'
253  << this->internalField().name() << " :"
254  << " heat transfer rate:" << Q
255  << " walltemperature "
256  << " min:" << gMin(Tp)
257  << " max:" << gMax(Tp)
258  << " avg:" << gAverage(Tp)
259  << endl;
260  }
261 
262  // Restore tag
263  UPstream::msgType() = oldTag;
264 }
265 
266 
268 (
269  Ostream& os
270 ) const
271 {
273  writeEntry(os, "Tnbr", TnbrName_);
274  writeEntry(os, "qrNbr", qrNbrName_);
275  writeEntry(os, "qr", qrName_);
276  writeEntry(os, "thicknessLayers", thicknessLayers_);
277  writeEntry(os, "kappaLayers", kappaLayers_);
278 
280 }
281 
282 
283 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
284 
286 (
289 );
290 
291 
292 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
293 
294 } // End namespace compressible
295 } // End namespace Foam
296 
297 
298 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:643
#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
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:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
Type gMin(const FieldField< Field, Type > &f)
const polyMesh & sampleMesh() const
Get the region mesh.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
Mixed boundary condition for temperature and radiation heat transfer to be used for in multiregion ca...
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
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
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.
turbulentTemperatureRadCoupledMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
faceListList boundary(nPatches)
virtual label size() const
Return size.
Definition: fvPatch.H:156
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
Type gMax(const FieldField< Field, Type > &f)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
const GeometricField::Patch & lookupPatchField(const word &name) const
Lookup and return the patchField of the named field from the.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
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,.
tmp< scalarField > kappa(const fvPatchScalarField &Tp) const
Given patch temperature calculate corresponding K field.
Type gAverage(const FieldField< Field, Type > &f)
const scalarField & deltaCoeffs() const
Return the face - cell distance coeffient.
Definition: fvPatch.C:170
label index() const
Return the index of this patch in the boundaryMesh.
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:275
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
virtual const word & name() const
Return name.
Definition: fvPatch.H:144
const polyPatch & samplePolyPatch() const
Get the patch on the region.
void write(Ostream &) const
Write.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
makePatchTypeField(fvPatchScalarField, thermalBaffleFvPatchScalarField)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844