turbulentTemperatureCoupledBaffleMixedFvPatchScalarField.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  thicknessLayers_(0),
52  kappaLayers_(0),
53  contactRes_(0)
54 {
55  this->refValue() = 0.0;
56  this->refGrad() = 0.0;
57  this->valueFraction() = 1.0;
58 }
59 
60 
63 (
65  const fvPatch& p,
67  const fvPatchFieldMapper& mapper
68 )
69 :
70  mixedFvPatchScalarField(ptf, p, iF, mapper),
71  temperatureCoupledBase(patch(), ptf),
72  TnbrName_(ptf.TnbrName_),
73  thicknessLayers_(ptf.thicknessLayers_),
74  kappaLayers_(ptf.kappaLayers_),
75  contactRes_(ptf.contactRes_)
76 {}
77 
78 
81 (
82  const fvPatch& p,
84  const dictionary& dict
85 )
86 :
87  mixedFvPatchScalarField(p, iF),
88  temperatureCoupledBase(patch(), dict),
89  TnbrName_(dict.lookup("Tnbr")),
90  thicknessLayers_(0),
91  kappaLayers_(0),
92  contactRes_(0.0)
93 {
94  if (!isA<mappedPatchBase>(this->patch().patch()))
95  {
97  (
98  "turbulentTemperatureCoupledBaffleMixedFvPatchScalarField::"
99  "turbulentTemperatureCoupledBaffleMixedFvPatchScalarField\n"
100  "(\n"
101  " const fvPatch& p,\n"
102  " const DimensionedField<scalar, volMesh>& iF,\n"
103  " const dictionary& dict\n"
104  ")\n"
105  ) << "\n patch type '" << p.type()
106  << "' not type '" << mappedPatchBase::typeName << "'"
107  << "\n for patch " << p.name()
108  << " of field " << dimensionedInternalField().name()
109  << " in file " << dimensionedInternalField().objectPath()
110  << exit(FatalError);
111  }
112 
113  if (dict.found("thicknessLayers"))
114  {
115  dict.lookup("thicknessLayers") >> thicknessLayers_;
116  dict.lookup("kappaLayers") >> kappaLayers_;
117 
118  if (thicknessLayers_.size() > 0)
119  {
120  // Calculate effective thermal resistance by harmonic averaging
121  forAll (thicknessLayers_, iLayer)
122  {
123  contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
124  }
125  contactRes_ = 1.0/contactRes_;
126  }
127  }
128 
129  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
130 
131  if (dict.found("refValue"))
132  {
133  // Full restart
134  refValue() = scalarField("refValue", dict, p.size());
135  refGrad() = scalarField("refGradient", dict, p.size());
136  valueFraction() = scalarField("valueFraction", dict, p.size());
137  }
138  else
139  {
140  // Start from user entered data. Assume fixedValue.
141  refValue() = *this;
142  refGrad() = 0.0;
143  valueFraction() = 1.0;
144  }
145 }
146 
147 
150 (
153 )
154 :
155  mixedFvPatchScalarField(wtcsf, iF),
156  temperatureCoupledBase(patch(), wtcsf),
157  TnbrName_(wtcsf.TnbrName_),
158  thicknessLayers_(wtcsf.thicknessLayers_),
159  kappaLayers_(wtcsf.kappaLayers_),
160  contactRes_(wtcsf.contactRes_)
161 {}
162 
163 
164 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
165 
167 {
168  if (updated())
169  {
170  return;
171  }
172 
173  // Since we're inside initEvaluate/evaluate there might be processor
174  // comms underway. Change the tag we use.
175  int oldTag = UPstream::msgType();
176  UPstream::msgType() = oldTag+1;
177 
178  // Get the coupling information from the mappedPatchBase
179  const mappedPatchBase& mpp =
180  refCast<const mappedPatchBase>(patch().patch());
181  const polyMesh& nbrMesh = mpp.sampleMesh();
182  const label samplePatchI = mpp.samplePolyPatch().index();
183  const fvPatch& nbrPatch =
184  refCast<const fvMesh>(nbrMesh).boundary()[samplePatchI];
185 
186  // Calculate the temperature by harmonic averaging
187  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
188 
190  refCast
191  <
193  >
194  (
195  nbrPatch.lookupPatchField<volScalarField, scalar>
196  (
197  TnbrName_
198  )
199  );
200 
201  // Swap to obtain full local values of neighbour internal field
202  tmp<scalarField> nbrIntFld(new scalarField(nbrField.size(), 0.0));
203  tmp<scalarField> nbrKDelta(new scalarField(nbrField.size(), 0.0));
204 
205  if (contactRes_ == 0.0)
206  {
207  nbrIntFld() = nbrField.patchInternalField();
208  nbrKDelta() = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
209  }
210  else
211  {
212  nbrIntFld() = nbrField;
213  nbrKDelta() = contactRes_;
214  }
215 
216  mpp.distribute(nbrIntFld());
217  mpp.distribute(nbrKDelta());
218 
219  tmp<scalarField> myKDelta = kappa(*this)*patch().deltaCoeffs();
220 
221 
222  // Both sides agree on
223  // - temperature : (myKDelta*fld + nbrKDelta*nbrFld)/(myKDelta+nbrKDelta)
224  // - gradient : (temperature-fld)*delta
225  // We've got a degree of freedom in how to implement this in a mixed bc.
226  // (what gradient, what fixedValue and mixing coefficient)
227  // Two reasonable choices:
228  // 1. specify above temperature on one side (preferentially the high side)
229  // and above gradient on the other. So this will switch between pure
230  // fixedvalue and pure fixedgradient
231  // 2. specify gradient and temperature such that the equations are the
232  // same on both sides. This leads to the choice of
233  // - refGradient = zero gradient
234  // - refValue = neighbour value
235  // - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())
236 
237  this->refValue() = nbrIntFld();
238  this->refGrad() = 0.0;
239  this->valueFraction() = nbrKDelta()/(nbrKDelta() + myKDelta());
240 
241  mixedFvPatchScalarField::updateCoeffs();
242 
243  if (debug)
244  {
245  scalar Q = gSum(kappa(*this)*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(*this)
256  << " max:" << gMax(*this)
257  << " avg:" << gAverage(*this)
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_
273  << token::END_STATEMENT << nl;
274  thicknessLayers_.writeEntry("thicknessLayers", os);
275  kappaLayers_.writeEntry("kappaLayers", os);
276 
278 }
279 
280 
281 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
282 
284 (
287 );
288 
289 
290 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
291 
292 } // End namespace compressible
293 } // End namespace Foam
294 
295 
296 // ************************************************************************* //
Mixed boundary condition for temperature, to be used for heat-transfer on back-to-back baffles...
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)
turbulentTemperatureCoupledBaffleMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
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.
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.
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
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
#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)
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)
A class for managing temporary objects.
Definition: PtrList.H:118