turbulentTemperatureCoupledBaffleMixedFvPatchScalarField.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  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 (
64  const fvPatch& p,
66  const dictionary& dict
67 )
68 :
69  mixedFvPatchScalarField(p, iF),
70  temperatureCoupledBase(patch(), dict),
71  TnbrName_(dict.lookup("Tnbr")),
72  thicknessLayers_(0),
73  kappaLayers_(0),
74  contactRes_(0.0)
75 {
76  if (!isA<mappedPatchBase>(this->patch().patch()))
77  {
79  << "' not type '" << mappedPatchBase::typeName << "'"
80  << "\n for patch " << p.name()
81  << " of field " << internalField().name()
82  << " in file " << internalField().objectPath()
83  << exit(FatalError);
84  }
85 
86  if (dict.found("thicknessLayers"))
87  {
88  dict.lookup("thicknessLayers") >> thicknessLayers_;
89  dict.lookup("kappaLayers") >> kappaLayers_;
90 
91  if (thicknessLayers_.size() > 0)
92  {
93  // Calculate effective thermal resistance by harmonic averaging
94  forAll(thicknessLayers_, iLayer)
95  {
96  contactRes_ += thicknessLayers_[iLayer]/kappaLayers_[iLayer];
97  }
98  contactRes_ = 1.0/contactRes_;
99  }
100  }
101 
102  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
103 
104  if (dict.found("refValue"))
105  {
106  // Full restart
107  refValue() = scalarField("refValue", dict, p.size());
108  refGrad() = scalarField("refGradient", dict, p.size());
109  valueFraction() = scalarField("valueFraction", dict, p.size());
110  }
111  else
112  {
113  // Start from user entered data. Assume fixedValue.
114  refValue() = *this;
115  refGrad() = 0.0;
116  valueFraction() = 1.0;
117  }
118 }
119 
120 
123 (
125  const fvPatch& p,
127  const fvPatchFieldMapper& mapper
128 )
129 :
130  mixedFvPatchScalarField(ptf, p, iF, mapper),
131  temperatureCoupledBase(patch(), ptf),
132  TnbrName_(ptf.TnbrName_),
133  thicknessLayers_(ptf.thicknessLayers_),
134  kappaLayers_(ptf.kappaLayers_),
135  contactRes_(ptf.contactRes_)
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 
182  const fvPatchScalarField& nbrTp =
183  nbrPatch.lookupPatchField<volScalarField, scalar>(TnbrName_);
184 
185  if (!isA<thisType>(nbrTp))
186  {
188  << "Patch field for " << internalField().name() << " on "
189  << patch().name() << " is of type " << thisType::typeName
190  << endl << "The neighbouring patch field " << TnbrName_ << " on "
191  << nbrPatch.name() << " is required to be the same, but is "
192  << "currently of type " << nbrTp.type() << exit(FatalError);
193  }
194 
195  const thisType& nbrField = refCast<const thisType>(nbrTp);
196 
197  // Swap to obtain full local values of neighbour internal field
198  tmp<scalarField> nbrIntFld(new scalarField(nbrField.size(), 0.0));
199  tmp<scalarField> nbrKDelta(new scalarField(nbrField.size(), 0.0));
200 
201  if (contactRes_ == 0.0)
202  {
203  nbrIntFld.ref() = nbrField.patchInternalField();
204  nbrKDelta.ref() = nbrField.kappa(nbrField)*nbrPatch.deltaCoeffs();
205  }
206  else
207  {
208  nbrIntFld.ref() = nbrField;
209  nbrKDelta.ref() = contactRes_;
210  }
211 
212  mpp.distribute(nbrIntFld.ref());
213  mpp.distribute(nbrKDelta.ref());
214 
215  tmp<scalarField> myKDelta = kappa(*this)*patch().deltaCoeffs();
216 
217 
218  // Both sides agree on
219  // - temperature : (myKDelta*fld + nbrKDelta*nbrFld)/(myKDelta+nbrKDelta)
220  // - gradient : (temperature-fld)*delta
221  // We've got a degree of freedom in how to implement this in a mixed bc.
222  // (what gradient, what fixedValue and mixing coefficient)
223  // Two reasonable choices:
224  // 1. specify above temperature on one side (preferentially the high side)
225  // and above gradient on the other. So this will switch between pure
226  // fixedvalue and pure fixedgradient
227  // 2. specify gradient and temperature such that the equations are the
228  // same on both sides. This leads to the choice of
229  // - refGradient = zero gradient
230  // - refValue = neighbour value
231  // - mixFraction = nbrKDelta / (nbrKDelta + myKDelta())
232 
233  this->refValue() = nbrIntFld();
234  this->refGrad() = 0.0;
235  this->valueFraction() = nbrKDelta()/(nbrKDelta() + myKDelta());
236 
237  mixedFvPatchScalarField::updateCoeffs();
238 
239  if (debug)
240  {
241  scalar Q = gSum(kappa(*this)*patch().magSf()*snGrad());
242 
243  Info<< patch().boundaryMesh().mesh().name() << ':'
244  << patch().name() << ':'
245  << this->internalField().name() << " <- "
246  << nbrMesh.name() << ':'
247  << nbrPatch.name() << ':'
248  << this->internalField().name() << " :"
249  << " heat transfer rate:" << Q
250  << " walltemperature "
251  << " min:" << gMin(*this)
252  << " max:" << gMax(*this)
253  << " avg:" << gAverage(*this)
254  << endl;
255  }
256 
257  // Restore tag
258  UPstream::msgType() = oldTag;
259 }
260 
261 
263 (
264  Ostream& os
265 ) const
266 {
268  writeEntry(os, "Tnbr", TnbrName_);
269  writeEntry(os, "thicknessLayers", thicknessLayers_);
270  writeEntry(os, "kappaLayers", kappaLayers_);
271 
273 }
274 
275 
276 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
277 
279 (
282 );
283 
284 
285 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
286 
287 } // End namespace compressible
288 } // End namespace Foam
289 
290 
291 // ************************************************************************* //
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: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
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
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.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
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
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)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
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.
Common functions used in temperature coupled boundaries.
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
A class for managing temporary objects.
Definition: PtrList.H:53
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