thermalBaffle.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "thermalBaffle.H"
30 #include "fvm.H"
31 #include "fvcDiv.H"
34 #include "fvMatrices.H"
36 
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41 namespace regionModels
42 {
43 namespace thermalBaffleModels
44 {
45 
46 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
47 
48 defineTypeNameAndDebug(thermalBaffle, 0);
49 
50 addToRunTimeSelectionTable(thermalBaffleModel, thermalBaffle, mesh);
51 addToRunTimeSelectionTable(thermalBaffleModel, thermalBaffle, dictionary);
52 
53 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
54 
56 {
57  this->solution().readEntry("nNonOrthCorr", nNonOrthCorr_);
58  return regionModel1D::read();
59 }
60 
61 
63 {
64  this->solution().readEntry("nNonOrthCorr", nNonOrthCorr_);
65  return regionModel1D::read(dict);
66 }
67 
68 
70 {
72 
73  const polyBoundaryMesh& rbm = regionMesh().boundaryMesh();
74 
76  (
77  new volScalarField
78  (
79  IOobject
80  (
81  "tQ",
82  regionMesh().time().timeName(),
83  regionMesh(),
86  false
87  ),
88  regionMesh(),
90  )
91  );
92 
93  volScalarField& Q = tQ.ref();
94 
95  volScalarField rho("rho", thermo_->rho());
96  volScalarField alpha("alpha", thermo_->alpha());
97 
98 
99  // If region is one-dimension variable thickness
100  if (oneD_ && !constantThickness_)
101  {
102  // Scale K and rhoCp and fill Q in the internal baffle region.
103  const label patchi = intCoupledPatchIDs_[0];
104  const polyPatch& ppCoupled = rbm[patchi];
105 
106  forAll(ppCoupled, localFacei)
107  {
108  const labelList& cells = boundaryFaceCells_[localFacei];
109  forAll(cells, i)
110  {
111  const label cellId = cells[i];
112 
113  Q[cellId] =
114  qs_.boundaryField()[patchi][localFacei]
115  /thickness_[localFacei];
116 
117  rho[cellId] *= delta_.value()/thickness_[localFacei];
118 
119  alpha[cellId] *= delta_.value()/thickness_[localFacei];
120  }
121  }
122  }
123  else
124  {
125  Q = Q_;
126  }
127 
128  fvScalarMatrix hEqn
129  (
130  fvm::ddt(rho, h_)
132  ==
133  Q
134  );
135 
136  if (moveMesh_)
137  {
138  surfaceScalarField phiMesh
139  (
141  );
142 
143  hEqn -= fvc::div(phiMesh);
144  }
145 
146  hEqn.relax();
147  hEqn.solve();
148 
149  thermo_->correct();
150 
151  Info<< "T min/max = " << min(thermo_->T()) << ", "
152  << max(thermo_->T()) << endl;
153 }
154 
155 
156 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
157 
158 thermalBaffle::thermalBaffle
159 (
160  const word& modelType,
161  const fvMesh& mesh,
162  const dictionary& dict
163 )
164 :
165  thermalBaffleModel(modelType, mesh, dict),
166  nNonOrthCorr_(solution().get<label>("nNonOrthCorr")),
167  thermo_(solidThermo::New(regionMesh(), dict)),
168  h_(thermo_->he()),
169  qs_
170  (
171  IOobject
172  (
173  "qs",
174  regionMesh().time().timeName(),
175  regionMesh(),
178  ),
179  regionMesh(),
181  ),
182  Q_
183  (
184  IOobject
185  (
186  "Q",
187  regionMesh().time().timeName(),
188  regionMesh(),
191  ),
192  regionMesh(),
194  ),
195  radiation_
196  (
198  (
199  dict.subDict("radiation"),
200  thermo_->T()
201  )
202  )
203 {
204  init();
205  thermo_->correct();
206 }
207 
208 
209 thermalBaffle::thermalBaffle
210 (
211  const word& modelType,
212  const fvMesh& mesh
213 )
214 :
215  thermalBaffleModel(modelType, mesh),
216  nNonOrthCorr_(solution().get<label>("nNonOrthCorr")),
217  thermo_(solidThermo::New(regionMesh())),
218  h_(thermo_->he()),
219  qs_
220  (
221  IOobject
222  (
223  "qs",
224  regionMesh().time().timeName(),
225  regionMesh(),
228  ),
229  regionMesh(),
231  ),
232  Q_
233  (
234  IOobject
235  (
236  "Q",
237  regionMesh().time().timeName(),
238  regionMesh(),
241  ),
242  regionMesh(),
244  ),
245  radiation_
246  (
248  (
249  thermo_->T()
250  )
251  )
252 {
253  init();
254  thermo_->correct();
255 }
256 
257 
258 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
259 
261 {}
262 
263 
264 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
265 
266 void thermalBaffle::init()
267 {
268  if (oneD_ && !constantThickness_)
269  {
270  label patchi = intCoupledPatchIDs_[0];
271  const label qsb = qs_.boundaryField()[patchi].size();
272 
273  if (qsb!= thickness_.size())
274  {
276  << "the boundary field of qs is "
277  << qsb << " and " << nl
278  << "the field 'thickness' is " << thickness_.size() << nl
279  << exit(FatalError);
280  }
281  }
282 }
283 
284 
286 {}
287 
288 
290 {
291  for (int nonOrth=0; nonOrth<=nNonOrthCorr_; nonOrth++)
292  {
293  solveEnergy();
294  }
295 }
296 
297 
299 {
300  return thermo_->Cp();
301 }
302 
303 
305 {
306  return radiation_->absorptionEmission().a();
307 }
308 
309 
311 {
312  return thermo_->rho();
313 }
314 
315 
317 {
318  return thermo_->kappa();
319 }
320 
321 
323 {
324  return thermo_->T();
325 }
326 
327 
329 {
330  return *thermo_;
331 }
332 
333 
335 {
336  const labelList& coupledPatches = intCoupledPatchIDs();
337 
338  forAll(coupledPatches, i)
339  {
340  const label patchi = coupledPatches[i];
341  const fvPatchScalarField& ph = h_.boundaryField()[patchi];
342  const word patchName = regionMesh().boundary()[patchi].name();
343  Info<< indent << "Q : " << patchName << indent <<
344  gSum
345  (
346  mag(regionMesh().Sf().boundaryField()[patchi])
347  * ph.snGrad()
348  * thermo_->alpha().boundaryField()[patchi]
349  ) << endl;
350  }
351 }
352 
353 
354 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
355 
356 } // End namespace thermalBaffleModels
357 } // End namespace regionModels
358 } // End namespace Foam
359 
360 // ************************************************************************* //
Foam::fvPatchField< scalar >
Foam::regionModels::thermalBaffleModels::thermalBaffleModel::oneD_
bool oneD_
Is it one dimension.
Definition: thermalBaffleModel.H:89
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::regionModels::regionModel1D::read
virtual bool read()
Read control parameters from dictionary.
Definition: regionModel1D.C:170
Foam::regionModels::thermalBaffleModels::thermalBaffle::T
virtual const volScalarField & T() const
Return temperature [K].
Definition: thermalBaffle.C:322
Foam::fvPatchField::snGrad
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:225
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::IOobject::AUTO_WRITE
Definition: IOobject.H:194
Foam::solution
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:55
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::regionModels::thermalBaffleModels::thermalBaffle::kappa
virtual const volScalarField & kappa() const
Return thermal conductivity [W/m/K].
Definition: thermalBaffle.C:316
Foam::solidThermo::New
static autoPtr< solidThermo > New(const fvMesh &, const word &phaseName=word::null)
Return a pointer to a new solidThermo created from.
Definition: solidThermo.C:82
Foam::polyBoundaryMesh
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
Definition: polyBoundaryMesh.H:63
Foam::regionModels::thermalBaffleModels::thermalBaffleModel::delta_
dimensionedScalar delta_
Baffle mesh thickness.
Definition: thermalBaffleModel.H:86
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::dimEnergy
const dimensionSet dimEnergy
Foam::regionModels::thermalBaffleModels::thermalBaffleModel::constantThickness_
bool constantThickness_
Is thickness constant.
Definition: thermalBaffleModel.H:92
Foam::regionModels::thermalBaffleModels::thermalBaffle::qs_
volScalarField qs_
Surface energy source / [J/m2/s].
Definition: thermalBaffle.H:98
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::regionModels::thermalBaffleModels::thermalBaffle::h_
volScalarField & h_
Enthalpy/internal energy.
Definition: thermalBaffle.H:92
thermalBaffle.H
fvcDiv.H
Calculate the divergence of the given field.
Foam::regionModels::thermalBaffleModels::thermalBaffle::~thermalBaffle
virtual ~thermalBaffle()
Destructor.
Definition: thermalBaffle.C:260
Foam::regionModels::regionModel::time
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:40
Foam::fvc::div
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:49
Foam::regionModels::thermalBaffleModels::thermalBaffle::solveEnergy
void solveEnergy()
Solve energy equation.
Definition: thermalBaffle.C:69
Foam::polyMesh::boundaryMesh
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:444
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::regionModels::thermalBaffleModels::thermalBaffle::Q_
volScalarField Q_
Volumetric energy source / [J/m3/s].
Definition: thermalBaffle.H:101
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::regionModels::thermalBaffleModels::thermalBaffle::thermo
virtual const solidThermo & thermo() const
Return const reference to the solidThermo.
Definition: thermalBaffle.C:328
Foam::regionModels::thermalBaffleModels::thermalBaffleModel::thickness_
scalarField thickness_
Baffle physical thickness.
Definition: thermalBaffleModel.H:83
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::regionModels::regionModel1D::boundaryFaceCells_
labelListList boundaryFaceCells_
Global cell IDs.
Definition: regionModel1D.H:83
fvMatrices.H
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Foam::solidThermo
Fundamental solid thermodynamic properties.
Definition: solidThermo.H:52
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::regionModels::thermalBaffleModels::thermalBaffleModel
Definition: thermalBaffleModel.H:60
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::radiation::radiationModel::New
static autoPtr< radiationModel > New(const volScalarField &T)
Return a reference to the selected radiation model.
Definition: radiationModelNew.C:36
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::fvm::laplacian
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:48
Foam::fvMatrix::relax
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:1183
Foam::regionModels::regionModel::regionMesh
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:64
init
mesh init(true)
Foam::IOobject::READ_IF_PRESENT
Definition: IOobject.H:187
Foam::regionModels::thermalBaffleModels::addToRunTimeSelectionTable
addToRunTimeSelectionTable(thermalBaffleModel, noThermo, mesh)
Foam::dictionary::readEntry
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
Definition: dictionaryTemplates.C:302
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::regionModels::regionModel::solution
const dictionary & solution() const
Return the solution dictionary.
Definition: regionModelI.H:99
timeName
word timeName
Definition: getTimeIndex.H:3
dict
dictionary dict
Definition: searchingEngine.H:14
fvm.H
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::fvc::interpolate
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
Foam::regionModels::thermalBaffleModels::thermalBaffle::evolveRegion
virtual void evolveRegion()
Evolve the thermal baffle.
Definition: thermalBaffle.C:289
Foam::indent
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:339
Foam::regionModels::thermalBaffleModels::thermalBaffle::radiation_
autoPtr< radiation::radiationModel > radiation_
Pointer to radiation model.
Definition: thermalBaffle.H:107
cellId
label cellId
Definition: interrogateWallPatches.H:67
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::regionModels::thermalBaffleModels::thermalBaffle::kappaRad
virtual const volScalarField & kappaRad() const
Return solid absortivity [1/m].
Definition: thermalBaffle.C:304
Foam::regionModels::thermalBaffleModels::thermalBaffle::preEvolveRegion
virtual void preEvolveRegion()
Pre-evolve thermal baffle.
Definition: thermalBaffle.C:285
Foam::fvMesh::boundary
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:685
Foam::regionModels::thermalBaffleModels::defineTypeNameAndDebug
defineTypeNameAndDebug(noThermo, 0)
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::fvm::ddt
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:48
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::regionModels::thermalBaffleModels::thermalBaffle::rho
virtual const volScalarField & rho() const
Return density [Kg/m3].
Definition: thermalBaffle.C:310
Foam::regionModels::thermalBaffleModels::thermalBaffle::info
virtual void info()
Provide some feedback.
Definition: thermalBaffle.C:334
Foam::regionModels::thermalBaffleModels::thermalBaffle::Cp
virtual const tmp< volScalarField > Cp() const
Return the film specific heat capacity [J/kg/K].
Definition: thermalBaffle.C:298
absorptionEmissionModel.H
Foam::List< label >
Foam::regionModels::thermalBaffleModels::thermalBaffle::read
virtual bool read()
Read control parameters IO dictionary.
Definition: thermalBaffle.C:55
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::fvMatrix
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvPatchField.H:68
Foam::fvMatrix::solve
SolverPerformance< Type > solve(const dictionary &)
Solve returning the solution statistics.
Definition: fvMatrixSolve.C:319
Foam::regionModels::regionModel1D::moveMesh_
Switch moveMesh_
Flag to allow mesh movement.
Definition: regionModel1D.H:98
Foam::regionModels::regionModel::intCoupledPatchIDs_
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
Definition: regionModel.H:114
cells
const cellShapeList & cells
Definition: gmvOutputHeader.H:3
Foam::regionModels::thermalBaffleModels::thermalBaffle::thermo_
autoPtr< solidThermo > thermo_
Solid thermo.
Definition: thermalBaffle.H:89
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
zeroGradientFvPatchFields.H
Foam::regionModels::regionModel::intCoupledPatchIDs
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
Definition: regionModelI.H:175
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::regionModels::thermalBaffleModels::thermalBaffle::nNonOrthCorr_
label nNonOrthCorr_
Number of non orthogonal correctors.
Definition: thermalBaffle.H:83