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-2016 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  << "' not type '" << mappedPatchBase::typeName << "'"
98  << "\n for patch " << p.name()
99  << " of field " << internalField().name()
100  << " in file " << internalField().objectPath()
101  << exit(FatalError);
102  }
103 
104  if (dict.found("thicknessLayers"))
105  {
106  dict.lookup("thicknessLayers") >> thicknessLayers_;
107  dict.lookup("kappaLayers") >> kappaLayers_;
108 
109  if (thicknessLayers_.size() > 0)
110  {
111  // Calculate effective thermal resistance by harmonic averaging
112  forAll(thicknessLayers_, iLayer)
113  {
114  contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
115  }
116  contactRes_ = 1.0/contactRes_;
117  }
118  }
119 
120  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
121 
122  if (dict.found("refValue"))
123  {
124  // Full restart
125  refValue() = scalarField("refValue", dict, p.size());
126  refGrad() = scalarField("refGradient", dict, p.size());
127  valueFraction() = scalarField("valueFraction", dict, p.size());
128  }
129  else
130  {
131  // Start from user entered data. Assume fixedValue.
132  refValue() = *this;
133  refGrad() = 0.0;
134  valueFraction() = 1.0;
135  }
136 }
137 
138 
141 (
144 )
145 :
146  mixedFvPatchScalarField(wtcsf, iF),
147  temperatureCoupledBase(patch(), wtcsf),
148  TnbrName_(wtcsf.TnbrName_),
149  thicknessLayers_(wtcsf.thicknessLayers_),
150  kappaLayers_(wtcsf.kappaLayers_),
151  contactRes_(wtcsf.contactRes_)
152 {}
153 
154 
155 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
156 
158 {
159  if (updated())
160  {
161  return;
162  }
163 
164  // Since we're inside initEvaluate/evaluate there might be processor
165  // comms underway. Change the tag we use.
166  int oldTag = UPstream::msgType();
167  UPstream::msgType() = oldTag+1;
168 
169  // Get the coupling information from the mappedPatchBase
170  const mappedPatchBase& mpp =
171  refCast<const mappedPatchBase>(patch().patch());
172  const polyMesh& nbrMesh = mpp.sampleMesh();
173  const label samplePatchi = mpp.samplePolyPatch().index();
174  const fvPatch& nbrPatch =
175  refCast<const fvMesh>(nbrMesh).boundary()[samplePatchi];
176 
177  // Calculate the temperature by harmonic averaging
178  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
179 
181  refCast
182  <
184  >
185  (
186  nbrPatch.lookupPatchField<volScalarField, scalar>
187  (
188  TnbrName_
189  )
190  );
191 
192  // Swap to obtain full local values of neighbour internal field
193  tmp<scalarField> nbrIntFld(new scalarField(nbrField.size(), 0.0));
194  tmp<scalarField> nbrKDelta(new scalarField(nbrField.size(), 0.0));
195 
196  if (contactRes_ == 0.0)
197  {
198  nbrIntFld.ref() = nbrField.patchInternalField();
199  nbrKDelta.ref() = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
200  }
201  else
202  {
203  nbrIntFld.ref() = nbrField;
204  nbrKDelta.ref() = contactRes_;
205  }
206 
207  mpp.distribute(nbrIntFld.ref());
208  mpp.distribute(nbrKDelta.ref());
209 
210  tmp<scalarField> myKDelta = kappa(*this)*patch().deltaCoeffs();
211 
212 
213  // Both sides agree on
214  // - temperature : (myKDelta*fld + nbrKDelta*nbrFld)/(myKDelta+nbrKDelta)
215  // - gradient : (temperature-fld)*delta
216  // We've got a degree of freedom in how to implement this in a mixed bc.
217  // (what gradient, what fixedValue and mixing coefficient)
218  // Two reasonable choices:
219  // 1. specify above temperature on one side (preferentially the high side)
220  // and above gradient on the other. So this will switch between pure
221  // fixedvalue and pure fixedgradient
222  // 2. specify gradient and temperature such that the equations are the
223  // same on both sides. This leads to the choice of
224  // - refGradient = zero gradient
225  // - refValue = neighbour value
226  // - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())
227 
228  this->refValue() = nbrIntFld();
229  this->refGrad() = 0.0;
230  this->valueFraction() = nbrKDelta()/(nbrKDelta() + myKDelta());
231 
232  mixedFvPatchScalarField::updateCoeffs();
233 
234  if (debug)
235  {
236  scalar Q = gSum(kappa(*this)*patch().magSf()*snGrad());
237 
238  Info<< patch().boundaryMesh().mesh().name() << ':'
239  << patch().name() << ':'
240  << this->internalField().name() << " <- "
241  << nbrMesh.name() << ':'
242  << nbrPatch.name() << ':'
243  << this->internalField().name() << " :"
244  << " heat transfer rate:" << Q
245  << " walltemperature "
246  << " min:" << gMin(*this)
247  << " max:" << gMax(*this)
248  << " avg:" << gAverage(*this)
249  << endl;
250  }
251 
252  // Restore tag
253  UPstream::msgType() = oldTag;
254 }
255 
256 
258 (
259  Ostream& os
260 ) const
261 {
263  os.writeKeyword("Tnbr")<< TnbrName_
264  << token::END_STATEMENT << nl;
265  thicknessLayers_.writeEntry("thicknessLayers", os);
266  kappaLayers_.writeEntry("kappaLayers", os);
267 
269 }
270 
271 
272 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
273 
275 (
278 );
279 
280 
281 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
282 
283 } // End namespace compressible
284 } // End namespace Foam
285 
286 
287 // ************************************************************************* //
Mixed boundary condition for temperature, to be used for heat-transfer on back-to-back baffles...
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
const scalarField & deltaCoeffs() const
Return the face - cell distance coeffient.
Definition: fvPatch.C:164
turbulentTemperatureCoupledBaffleMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void write(Ostream &) const
Write.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
const word & name() const
Return name.
Definition: fvPatch.H:149
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
const polyMesh & sampleMesh() const
Get the region mesh.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:65
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:464
void writeEntry(Ostream &) const
Write the UList as a dictionary entry.
Definition: UListIO.C:35
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
Type gSum(const FieldField< Field, Type > &f)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
void distribute(List< Type > &lst) const
Wrapper around map/interpolate data distribution.
virtual label size() const
Return size.
Definition: fvPatch.H:161
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
faceListList boundary(nPatches)
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=NULL, const Type *=NULL) const
Lookup and return the patchField of the named field from the.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
Type gMax(const FieldField< Field, Type > &f)
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
Common functions used in temperature coupled boundaries.
const polyPatch & samplePolyPatch() const
Get the patch on the region.
Type gAverage(const FieldField< Field, Type > &f)
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:396
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
A class for managing temporary objects.
Definition: PtrList.H:54
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
tmp< scalarField > kappa(const scalarField &Tp) const
Given patch temperature calculate corresponding K field.
label index() const
Return the index of this patch in the boundaryMesh.
runTime 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:451
bool compressible
Definition: pEqn.H:30