thermalBaffle.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-2019 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 
26 #include "thermalBaffle.H"
27 
28 #include "fvm.H"
29 #include "fvcDiv.H"
32 #include "fvMatrices.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 namespace regionModels
40 {
41 namespace thermalBaffleModels
42 {
43 
44 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45 
46 defineTypeNameAndDebug(thermalBaffle, 0);
47 
48 addToRunTimeSelectionTable(thermalBaffleModel, thermalBaffle, mesh);
49 addToRunTimeSelectionTable(thermalBaffleModel, thermalBaffle, dictionary);
50 
51 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
52 
54 {
55  this->solution().lookup("nNonOrthCorr") >> nNonOrthCorr_;
56  return regionModel1D::read();
57 }
58 
59 
61 {
62  this->solution().lookup("nNonOrthCorr") >> nNonOrthCorr_;
63  return regionModel1D::read(dict);
64 }
65 
66 
68 {
70 
71  const polyBoundaryMesh& rbm = regionMesh().boundaryMesh();
72 
74  (
76  (
77  "tQ",
78  regionMesh(),
80  )
81  );
82 
83  volScalarField& Q = tQ.ref();
84 
85  volScalarField rho("rho", thermo_->rho());
86  volScalarField alpha("alpha", thermo_->alpha());
87 
88 
89  // If region is one-dimension variable thickness
90  if (oneD_ && !constantThickness_)
91  {
92  // Scale K and rhoCp and fill Q in the internal baffle region.
93  const label patchi = intCoupledPatchIDs_[0];
94  const polyPatch& ppCoupled = rbm[patchi];
95 
96  forAll(ppCoupled, localFacei)
97  {
98  const labelList& cells = boundaryFaceCells_[localFacei];
99  forAll(cells, i)
100  {
101  const label cellId = cells[i];
102 
103  Q[cellId] =
104  Qs_.boundaryField()[patchi][localFacei]
105  /thickness_[localFacei];
106 
107  rho[cellId] *= delta_.value()/thickness_[localFacei];
108 
109  alpha[cellId] *= delta_.value()/thickness_[localFacei];
110  }
111  }
112  }
113  else
114  {
115  Q = Q_;
116  }
117 
118  fvScalarMatrix hEqn
119  (
120  fvm::ddt(rho, h_)
122  ==
123  Q
124  );
125 
126  if (moveMesh_)
127  {
128  surfaceScalarField phiMesh
129  (
131  );
132 
133  hEqn -= fvc::div(phiMesh);
134  }
135 
136  hEqn.relax();
137  hEqn.solve();
138 
139  thermo_->correct();
140 
141  Info<< "T min/max = " << min(thermo_->T()) << ", "
142  << max(thermo_->T()) << endl;
143 }
144 
145 
146 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
147 
149 (
150  const word& modelType,
151  const fvMesh& mesh,
152  const dictionary& dict
153 )
154 :
155  thermalBaffleModel(modelType, mesh, dict),
156  nNonOrthCorr_(solution().lookup<label>("nNonOrthCorr")),
158  h_(thermo_->he()),
159  Qs_
160  (
161  IOobject
162  (
163  "Qs",
164  regionMesh().time().timeName(),
165  regionMesh(),
168  ),
169  regionMesh(),
171  ),
172  Q_
173  (
174  IOobject
175  (
176  "Q",
177  regionMesh().time().timeName(),
178  regionMesh(),
181  ),
182  regionMesh(),
184  ),
185  radiation_
186  (
188  (
189  dict.subDict("radiation"),
190  thermo_->T()
191  )
192  )
193 {
194  init();
195  thermo_->correct();
196 }
197 
198 
200 (
201  const word& modelType,
202  const fvMesh& mesh
203 )
204 :
205  thermalBaffleModel(modelType, mesh),
206  nNonOrthCorr_(solution().lookup<label>("nNonOrthCorr")),
208  h_(thermo_->he()),
209  Qs_
210  (
211  IOobject
212  (
213  "Qs",
214  regionMesh().time().timeName(),
215  regionMesh(),
218  ),
219  regionMesh(),
221  ),
222  Q_
223  (
224  IOobject
225  (
226  "Q",
227  regionMesh().time().timeName(),
228  regionMesh(),
231  ),
232  regionMesh(),
234  ),
235  radiation_
236  (
238  (
239  thermo_->T()
240  )
241  )
242 {
243  init();
244  thermo_->correct();
245 }
246 
247 
248 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
249 
251 {}
252 
253 
254 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
255 
256 void thermalBaffle::init()
257 {
258  if (oneD_ && !constantThickness_)
259  {
261  const label Qsb = Qs_.boundaryField()[patchi].size();
262 
263  if (Qsb!= thickness_.size())
264  {
266  << "the boundary field of Qs is "
267  << Qsb << " and " << nl
268  << "the field 'thickness' is " << thickness_.size() << nl
269  << exit(FatalError);
270  }
271  }
272 }
273 
274 
276 {}
277 
278 
280 {
281  for (int nonOrth=0; nonOrth<=nNonOrthCorr_; nonOrth++)
282  {
283  solveEnergy();
284  }
285 }
286 
287 
289 {
290  return thermo_->Cp();
291 }
292 
293 
295 {
296  return radiation_->absorptionEmission().a();
297 }
298 
299 
301 {
302  return thermo_->rho();
303 }
304 
305 
307 {
308  return thermo_->kappa();
309 }
310 
311 
313 {
314  return thermo_->T();
315 }
316 
317 
319 {
320  return thermo_;
321 }
322 
323 
325 {
326  const labelList& coupledPatches = intCoupledPatchIDs();
327 
328  forAll(coupledPatches, i)
329  {
330  const label patchi = coupledPatches[i];
331  const fvPatchScalarField& ph = h_.boundaryField()[patchi];
332  const word patchName = regionMesh().boundary()[patchi].name();
333  Info<< indent << "Q : " << patchName << indent <<
334  gSum
335  (
336  mag(regionMesh().Sf().boundaryField()[patchi])
337  * ph.snGrad()
338  * thermo_->alpha().boundaryField()[patchi]
339  ) << endl;
340  }
341 }
342 
343 
344 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
345 
346 } // end namespace thermalBaffleModels
347 } // end namespace regionModels
348 } // end namespace Foam
349 
350 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
autoPtr< solidThermo > thermo_
Solid thermo.
Definition: thermalBaffle.H:79
thermalBaffle(const word &modelType, const fvMesh &mesh)
Construct from components.
volScalarField Qs_
Surface energy source / [J/m2/s].
Definition: thermalBaffle.H:88
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
thermalBaffleModel(const fvMesh &mesh)
Construct null from mesh.
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 dimensionSet dimArea
autoPtr< radiationModel > radiation_
Pointer to radiation model.
Definition: thermalBaffle.H:97
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:221
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:166
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
label nNonOrthCorr_
Number of non orthogonal correctors.
Definition: thermalBaffle.H:73
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
virtual const solidThermo & thermo() const
Return const reference to the solidThermo.
static autoPtr< solidThermo > New(const fvMesh &, const word &phaseName=word::null)
Standard selection based on fvMesh.
Definition: solidThermo.C:113
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const Time & time() const
Return the reference to the time database.
Definition: regionModelI.H:37
static autoPtr< radiationModel > New(const volScalarField &T)
Return a reference to the selected radiation model.
Macros for easy insertion into run-time selection tables.
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:982
volScalarField & h_
Enthalpy/internal energy.
Definition: thermalBaffle.H:82
Switch moveMesh_
Flag to allow mesh movement.
Definition: regionModel1D.H:96
const dimensionSet dimTime
virtual void evolveRegion()
Evolve the thermal baffle.
Type gSum(const FieldField< Field, Type > &f)
const cellShapeList & cells
virtual const volScalarField & kappaRad() const
Return solid absorptivity [1/m].
virtual void preEvolveRegion()
Pre-evolve thermal baffle.
volScalarField Q_
Volumetric energy source / [J/m3/s].
Definition: thermalBaffle.H:91
addToRunTimeSelectionTable(thermalBaffleModel, noThermo, mesh)
A class for handling words, derived from string.
Definition: word.H:59
#define DebugInFunction
Report an information message using Foam::Info.
labelListList boundaryFaceCells_
Global cell IDs.
Definition: regionModel1D.H:81
const dictionary & solution() const
Return the solution dictionary.
Definition: regionModelI.H:97
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
const Type & value() const
Return const reference to value.
const fvMesh & regionMesh() const
Return the region mesh database.
Definition: regionModelI.H:55
phi
Definition: correctPhi.H:3
A special matrix type and solver, designed for finite volume solutions of scalar equations. Face addressing is used to make all matrix assembly and solution loops vectorise.
Definition: fvPatchField.H:72
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
virtual const tmp< volScalarField > Cp() const
Return the film specific heat capacity [J/kg/K].
Foam::polyBoundaryMesh.
Calculate the divergence of the given field.
virtual bool read()
Read control parameters IO dictionary.
Definition: thermalBaffle.C:53
SolverPerformance< Type > solve(const dictionary &)
Solve segregated or coupled returning the solution statistics.
Definition: fvMatrixSolve.C:58
static const char nl
Definition: Ostream.H:260
Base-class for solid thermodynamic properties.
Definition: solidThermo.H:53
void relax(const scalar alpha)
Relax matrix (for steady-state solution).
Definition: fvMatrix.C:521
const dimensionSet dimEnergy
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
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.
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
virtual bool read()
Read control parameters from dictionary.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
label cellId
virtual const volScalarField & T() const
Return temperature [K].
A special matrix type and solver, designed for finite volume solutions of scalar equations.
const dimensionSet dimVolume
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
labelList intCoupledPatchIDs_
List of patch IDs internally coupled with the primary region.
Definition: regionModel.H:106
A class for managing temporary objects.
Definition: PtrList.H:53
virtual const volScalarField & kappa() const
Return thermal conductivity [W/m/K].
const labelList & intCoupledPatchIDs() const
Return the list of patch IDs internally coupled with the.
Definition: regionModelI.H:173
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual void info()
Provide some feedback.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual const volScalarField & rho() const
Return density [Kg/m^3].
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
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540