turbulentTemperatureRadCoupledMixedFvPatchScalarField.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) 2011-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 
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(), "undefined", "undefined", "undefined-K"),
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 (
67  const fvPatch& p,
69  const fvPatchFieldMapper& mapper
70 )
71 :
72  mixedFvPatchScalarField(psf, p, iF, mapper),
73  temperatureCoupledBase(patch(), psf),
74  TnbrName_(psf.TnbrName_),
75  QrNbrName_(psf.QrNbrName_),
76  QrName_(psf.QrName_),
77  thicknessLayers_(psf.thicknessLayers_),
78  kappaLayers_(psf.kappaLayers_),
79  contactRes_(psf.contactRes_)
80 {}
81 
82 
85 (
86  const fvPatch& p,
88  const dictionary& dict
89 )
90 :
91  mixedFvPatchScalarField(p, iF),
92  temperatureCoupledBase(patch(), dict),
93  TnbrName_(dict.lookupOrDefault<word>("Tnbr", "T")),
94  QrNbrName_(dict.lookupOrDefault<word>("QrNbr", "none")),
95  QrName_(dict.lookupOrDefault<word>("Qr", "none")),
96  thicknessLayers_(0),
97  kappaLayers_(0),
98  contactRes_(0.0)
99 {
100  if (!isA<mappedPatchBase>(this->patch().patch()))
101  {
103  (
104  "turbulentTemperatureRadCoupledMixedFvPatchScalarField::"
105  "turbulentTemperatureRadCoupledMixedFvPatchScalarField\n"
106  "(\n"
107  " const fvPatch& p,\n"
108  " const DimensionedField<scalar, volMesh>& iF,\n"
109  " const dictionary& dict\n"
110  ")\n"
111  ) << "\n patch type '" << p.type()
112  << "' not type '" << mappedPatchBase::typeName << "'"
113  << "\n for patch " << p.name()
114  << " of field " << dimensionedInternalField().name()
115  << " in file " << dimensionedInternalField().objectPath()
116  << exit(FatalError);
117  }
118 
119  if (dict.found("thicknessLayers"))
120  {
121  dict.lookup("thicknessLayers") >> thicknessLayers_;
122  dict.lookup("kappaLayers") >> kappaLayers_;
123 
124  if (thicknessLayers_.size() > 0)
125  {
126  // Calculate effective thermal resistance by harmonic averaging
127  forAll (thicknessLayers_, iLayer)
128  {
129  contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
130  }
131  contactRes_ = 1.0/contactRes_;
132  }
133  }
134 
135  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
136 
137  if (dict.found("refValue"))
138  {
139  // Full restart
140  refValue() = scalarField("refValue", dict, p.size());
141  refGrad() = scalarField("refGradient", dict, p.size());
142  valueFraction() = scalarField("valueFraction", dict, p.size());
143  }
144  else
145  {
146  // Start from user entered data. Assume fixedValue.
147  refValue() = *this;
148  refGrad() = 0.0;
149  valueFraction() = 1.0;
150  }
151 }
152 
153 
156 (
159 )
160 :
161  mixedFvPatchScalarField(psf, iF),
162  temperatureCoupledBase(patch(), psf),
163  TnbrName_(psf.TnbrName_),
164  QrNbrName_(psf.QrNbrName_),
165  QrName_(psf.QrName_),
166  thicknessLayers_(psf.thicknessLayers_),
167  kappaLayers_(psf.kappaLayers_),
168  contactRes_(psf.contactRes_)
169 {}
170 
171 
172 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
173 
175 {
176  if (updated())
177  {
178  return;
179  }
180 
181  // Since we're inside initEvaluate/evaluate there might be processor
182  // comms underway. Change the tag we use.
183  int oldTag = UPstream::msgType();
184  UPstream::msgType() = oldTag+1;
185 
186  // Get the coupling information from the mappedPatchBase
187  const mappedPatchBase& mpp =
188  refCast<const mappedPatchBase>(patch().patch());
189  const polyMesh& nbrMesh = mpp.sampleMesh();
190  const label samplePatchI = mpp.samplePolyPatch().index();
191  const fvPatch& nbrPatch =
192  refCast<const fvMesh>(nbrMesh).boundary()[samplePatchI];
193 
194 
195  scalarField Tc(patchInternalField());
196  scalarField& Tp = *this;
197 
199  nbrField = refCast
201  (
202  nbrPatch.lookupPatchField<volScalarField, scalar>(TnbrName_)
203  );
204 
205  // Swap to obtain full local values of neighbour internal field
206  scalarField TcNbr(nbrField.patchInternalField());
207  mpp.distribute(TcNbr);
208 
209 
210  // Swap to obtain full local values of neighbour K*delta
211  scalarField KDeltaNbr;
212  if (contactRes_ == 0.0)
213  {
214  KDeltaNbr = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
215  }
216  else
217  {
218  KDeltaNbr.setSize(nbrField.size(), contactRes_);
219  }
220  mpp.distribute(KDeltaNbr);
221 
222  scalarField KDelta(kappa(Tp)*patch().deltaCoeffs());
223 
224  scalarField Qr(Tp.size(), 0.0);
225  if (QrName_ != "none")
226  {
227  Qr = patch().lookupPatchField<volScalarField, scalar>(QrName_);
228  }
229 
230  scalarField QrNbr(Tp.size(), 0.0);
231  if (QrNbrName_ != "none")
232  {
233  QrNbr = nbrPatch.lookupPatchField<volScalarField, scalar>(QrNbrName_);
234  mpp.distribute(QrNbr);
235  }
236 
237  valueFraction() = KDeltaNbr/(KDeltaNbr + KDelta);
238  refValue() = TcNbr;
239  refGrad() = (Qr + QrNbr)/kappa(Tp);
240 
241  mixedFvPatchScalarField::updateCoeffs();
242 
243  if (debug)
244  {
245  scalar Q = gSum(kappa(Tp)*patch().magSf()*snGrad());
246 
247  Info<< patch().boundaryMesh().mesh().name() << ':'
248  << patch().name() << ':'
249  << this->dimensionedInternalField().name() << " <- "
250  << nbrMesh.name() << ':'
251  << nbrPatch.name() << ':'
252  << this->dimensionedInternalField().name() << " :"
253  << " heat transfer rate:" << Q
254  << " walltemperature "
255  << " min:" << gMin(Tp)
256  << " max:" << gMax(Tp)
257  << " avg:" << gAverage(Tp)
258  << endl;
259  }
260 
261  // Restore tag
262  UPstream::msgType() = oldTag;
263 }
264 
265 
267 (
268  Ostream& os
269 ) const
270 {
272  os.writeKeyword("Tnbr")<< TnbrName_ << token::END_STATEMENT << nl;
273  os.writeKeyword("QrNbr")<< QrNbrName_ << token::END_STATEMENT << nl;
274  os.writeKeyword("Qr")<< QrName_ << token::END_STATEMENT << nl;
275  thicknessLayers_.writeEntry("thicknessLayers", os);
276  kappaLayers_.writeEntry("kappaLayers", os);
277 
279 }
280 
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
285 (
288 );
289 
290 
291 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
292 
293 } // End namespace compressible
294 } // End namespace Foam
295 
296 
297 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:387
faceListList boundary(nPatches)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
turbulentTemperatureRadCoupledMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const polyMesh & sampleMesh() const
Get the region mesh.
void writeEntry(Ostream &) const
Write the UList as a dictionary entry.
Definition: UListIO.C:35
const scalarField & deltaCoeffs() const
Return the face - cell distance coeffient.
Definition: fvPatch.C:164
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
label index() const
Return the index of this patch in the boundaryMesh.
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
void write(Ostream &) const
Write.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Foam::fvPatchFieldMapper.
messageStream Info
makePatchTypeField(fvPatchScalarField, thermalBaffleFvPatchScalarField)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
const word & name() const
Return name.
Definition: fvPatch.H:149
runTime write()
const GeometricField::PatchFieldType & lookupPatchField(const word &name, const GeometricField *=NULL, const Type *=NULL) const
Lookup and return the patchField of the named field from the.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Type gSum(const FieldField< Field, Type > &f)
static const char nl
Definition: Ostream.H:260
void setSize(const label)
Reset size of List.
Definition: List.C:318
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const polyPatch & samplePolyPatch() const
Get the patch on the region.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
#define forAll(list, i)
Definition: UList.H:421
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:451
Macros for easy insertion into run-time selection tables.
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
Type gMin(const FieldField< Field, Type > &f)
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
Mixed boundary condition for temperature and radiation heat transfer to be used for in multiregion ca...
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
rDeltaT dimensionedInternalField()
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
bool compressible
Definition: pEqn.H:40
error FatalError
Type gAverage(const FieldField< Field, Type > &f)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
Common functions for use in temperature coupled boundaries.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
virtual label size() const
Return size.
Definition: fvPatch.H:161
This function object calculates and outputs the second invariant of the velocity gradient tensor [1/s...
Definition: Q.H:66
Type gMax(const FieldField< Field, Type > &f)