thermalBaffleFvPatchScalarField.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 
27 #include "emptyPolyPatch.H"
28 #include "mappedWallPolyPatch.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace compressible
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
42 (
43  const fvPatch& p,
45 )
46 :
48  owner_(false),
49  baffle_(),
50  dict_(dictionary::null),
51  extrudeMeshPtr_()
52 {}
53 
54 
57 (
58  const fvPatch& p,
60  const dictionary& dict
61 )
62 :
64  owner_(false),
65  baffle_(),
66  dict_(dict),
67  extrudeMeshPtr_()
68 {
69 
70  const fvMesh& thisMesh = patch().boundaryMesh().mesh();
71 
73 
74  if (thisMesh.name() == polyMesh::defaultRegion)
75  {
76  const word regionName =
77  dict_.lookupOrDefault<word>("regionName", "none");
78 
79  const word baffleName("3DBaffle" + regionName);
80 
81  if
82  (
83  !thisMesh.time().foundObject<fvMesh>(regionName)
84  && regionName != "none"
85  )
86  {
87  if (extrudeMeshPtr_.empty())
88  {
89  createPatchMesh();
90  }
91 
92  baffle_.reset(baffle::New(thisMesh, dict).ptr());
93  owner_ = true;
94  baffle_->rename(baffleName);
95  }
96  }
97 }
98 
99 
102 (
104  const fvPatch& p,
106  const fvPatchFieldMapper& mapper
107 )
108 :
110  (
111  ptf,
112  p,
113  iF,
114  mapper
115  ),
116  owner_(ptf.owner_),
117  baffle_(),
118  dict_(ptf.dict_),
119  extrudeMeshPtr_()
120 {}
121 
122 
125 (
128 )
129 :
131  owner_(ptf.owner_),
132  baffle_(),
133  dict_(ptf.dict_),
134  extrudeMeshPtr_()
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
140 
142 (
143  const fvPatchFieldMapper& m
144 )
145 {
146  mixedFvPatchScalarField::autoMap(m);
147 }
148 
149 
151 (
152  const fvPatchScalarField& ptf,
153  const labelList& addr
154 )
155 {
156  mixedFvPatchScalarField::rmap(ptf, addr);
157 }
158 
159 
160 void thermalBaffleFvPatchScalarField::createPatchMesh()
161 {
162 
163  const fvMesh& thisMesh = patch().boundaryMesh().mesh();
164 
165  word regionName = dict_.lookup("regionName");
166 
167  List<polyPatch*> regionPatches(3);
168  List<word> patchNames(regionPatches.size());
169  List<word> patchTypes(regionPatches.size());
170  List<dictionary> dicts(regionPatches.size());
171 
172  patchNames[bottomPatchID] = word("bottom");
173  patchNames[sidePatchID] = word("side");
174  patchNames[topPatchID] = word("top");
175 
176  patchTypes[bottomPatchID] = mappedWallPolyPatch::typeName;
177  patchTypes[topPatchID] = mappedWallPolyPatch::typeName;
178 
179  if (readBool(dict_.lookup("columnCells")))
180  {
181  patchTypes[sidePatchID] = emptyPolyPatch::typeName;
182  }
183  else
184  {
185  patchTypes[sidePatchID] = polyPatch::typeName;
186  }
187 
188  const mappedPatchBase& mpp =
189  refCast<const mappedPatchBase>(patch().patch());
190 
191  const word coupleGroup(mpp.coupleGroup());
192 
193  wordList inGroups(1);
194  inGroups[0] = coupleGroup;
195 
196  dicts[bottomPatchID].add("coupleGroup", coupleGroup);
197  dicts[bottomPatchID].add("inGroups", inGroups);
198  dicts[bottomPatchID].add("sampleMode", mpp.sampleModeNames_[mpp.mode()]);
199 
200  const label sepPos = coupleGroup.find('_');
201 
202  const word coupleGroupSlave = coupleGroup(0, sepPos) + "_slave";
203 
204  inGroups[0] = coupleGroupSlave;
205  dicts[topPatchID].add("coupleGroup", coupleGroupSlave);
206  dicts[topPatchID].add("inGroups", inGroups);
207  dicts[topPatchID].add("sampleMode", mpp.sampleModeNames_[mpp.mode()]);
208 
209 
210  forAll(regionPatches, patchi)
211  {
212  dictionary& patchDict = dicts[patchi];
213  patchDict.set("nFaces", 0);
214  patchDict.set("startFace", 0);
215 
216  regionPatches[patchi] = polyPatch::New
217  (
219  patchNames[patchi],
220  dicts[patchi],
221  patchi,
222  thisMesh.boundaryMesh()
223  ).ptr();
224  }
225 
226  extrudeMeshPtr_.reset
227  (
228  new extrudePatchMesh
229  (
230  thisMesh,
231  patch(),
232  dict_,
233  regionName,
234  regionPatches
235  )
236  );
237 
238  if (extrudeMeshPtr_.empty())
239  {
241  << "Specified IOobject::MUST_READ_IF_MODIFIED but class"
242  << " patchMeshPtr not set."
243  << endl;
244  }
245 }
246 
247 
249 {
250  if (this->updated())
251  {
252  return;
253  }
254 
255  const fvMesh& thisMesh = patch().boundaryMesh().mesh();
256 
257  if (owner_ && thisMesh.name() == polyMesh::defaultRegion)
258  {
259  baffle_->evolve();
260  }
261 
263 }
264 
265 
267 {
269 
270  const fvMesh& thisMesh = patch().boundaryMesh().mesh();
271 
272  if (thisMesh.name() == polyMesh::defaultRegion && owner_)
273  {
274 
275  writeKeyword(os, "extrudeModel");
276  os << word(dict_.lookup("extrudeModel"))
277  << token::END_STATEMENT << nl;
278 
279  writeKeyword(os, "nLayers");
280  os << dict_.lookup<label>("nLayers")
281  << token::END_STATEMENT << nl;
282 
283  writeKeyword(os, "expansionRatio");
284  os << dict_.lookup<scalar>("expansionRatio")
285  << token::END_STATEMENT << nl;
286 
287  writeKeyword(os, "columnCells");
288  os << readBool(dict_.lookup("columnCells"))
289  << token::END_STATEMENT << nl;
290 
291  word extrudeModel(word(dict_.lookup("extrudeModel")) + "Coeffs");
292  writeKeyword(os, extrudeModel);
293  os << dict_.subDict(extrudeModel) << nl;
294 
295  word regionName = dict_.lookup("regionName");
296  writeKeyword(os, "regionName") << regionName
297  << token::END_STATEMENT << nl;
298 
299  writeKeyword(os, "thermoType");
300  os << dict_.subDict("thermoType") << nl;
301 
302  writeKeyword(os, "mixture");
303  os << dict_.subDict("mixture") << nl;
304 
305  writeKeyword(os, "radiation");
306  os << dict_.subDict("radiation") << nl;
307  }
308 }
309 
310 
311 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
312 
314 (
317 );
318 
319 
320 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
321 
322 } // End namespace compressible
323 } // End namespace Foam
324 
325 
326 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
#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
Ostream & writeKeyword(Foam::Ostream &os, const keyType &kw)
Write the keyword to the Ostream with the current level of indentation.
Definition: keyType.C:155
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
Foam::word regionName
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
wordList patchTypes(nPatches)
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:309
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
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Mixed boundary condition for temperature and radiation heat transfer to be used for in multiregion ca...
bool foundObject(const word &name) const
Is the named Type found?
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
Mesh at a patch created on the fly. The following entry should be used on the field boundary dictiona...
bool readBool(Istream &)
Definition: boolIO.C:60
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
Macros for easy insertion into run-time selection tables.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Top level extrusion model class.
Definition: extrudeModel.H:51
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:982
const sampleMode & mode() const
What to sample.
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
wordList patchNames(nPatches)
static const dictionary null
Null dictionary.
Definition: dictionary.H:242
const polyMesh & mesh() const
Return the mesh reference.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
static const char nl
Definition: Ostream.H:260
static const NamedEnum< sampleMode, 6 > sampleModeNames_
thermalBaffleFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
const word & name() const
Return reference to name.
Definition: fvMesh.H:253
const word & coupleGroup() const
PatchGroup (only if NEARESTPATCHFACE)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void set(entry *)
Assign a new entry, overwrite any existing entry.
Definition: dictionary.C:1271
static autoPtr< polyPatch > New(const word &patchType, const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm)
Return a pointer to a new patch created on freestore from.
Definition: polyPatchNew.C:32
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
This boundary condition provides a coupled temperature condition between multiple mesh regions...
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