variableHeatTransfer.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-2015 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 "variableHeatTransfer.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace fv
35 {
36  defineTypeNameAndDebug(variableHeatTransfer, 0);
38  (
39  option,
40  variableHeatTransfer,
41  dictionary
42  );
43 }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
50 (
51  const word& name,
52  const word& modelType,
53  const dictionary& dict,
54  const fvMesh& mesh
55 )
56 :
57  interRegionHeatTransferModel(name, modelType, dict, mesh),
58  UNbrName_(coeffs_.lookupOrDefault<word>("UNbrName", "U")),
59  a_(0),
60  b_(0),
61  c_(0),
62  ds_(0),
63  Pr_(0),
64  AoV_()
65 {
66  if (master_)
67  {
68  a_ = readScalar(coeffs_.lookup("a"));
69  b_ = readScalar(coeffs_.lookup("b"));
70  c_ = readScalar(coeffs_.lookup("c"));
71  ds_ = readScalar(coeffs_.lookup("ds"));
72  Pr_ = readScalar(coeffs_.lookup("Pr"));
73  AoV_.reset
74  (
75  new volScalarField
76  (
77  IOobject
78  (
79  "AoV",
80  mesh_.time().timeName(),
81  mesh_,
84  ),
85  mesh_
86  )
87  );
88  }
89 }
90 
91 
92 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
93 
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
101 {
102  const fvMesh& nbrMesh =
103  mesh_.time().lookupObject<fvMesh>(nbrRegionName());
104 
105  const compressible::turbulenceModel& nbrTurb =
107  (
109  );
110 
111  const fluidThermo& nbrThermo =
113 
114  const volVectorField& UNbr =
115  nbrMesh.lookupObject<volVectorField>(UNbrName_);
116 
117  const volScalarField ReNbr(mag(UNbr)*ds_*nbrThermo.rho()/nbrTurb.mut());
118 
119  const volScalarField NuNbr(a_*pow(ReNbr, b_)*pow(Pr_, c_));
120 
121  const scalarField htcNbr(NuNbr*nbrTurb.kappaEff()/ds_);
122 
123  const scalarField htcNbrMapped(interpolate(htcNbr));
124 
125  htc_.internalField() = htcNbrMapped*AoV_;
126 }
127 
128 
130 {
132  {
133  coeffs_.readIfPresent("UNbrName", UNbrName_);
134 
135  coeffs_.readIfPresent("a", a_);
136  coeffs_.readIfPresent("b", b_);
137  coeffs_.readIfPresent("c", c_);
138  coeffs_.readIfPresent("ds", ds_);
139  coeffs_.readIfPresent("Pr", Pr_);
140 
141  return true;
142  }
143  else
144  {
145  return false;
146  }
147 }
148 
149 
150 // ************************************************************************* //
defineTypeNameAndDebug(cellSetOption, 0)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual bool read(const dictionary &dict)
Read dictionary.
virtual void calculateHtc()
Calculate the heat transfer coefficient.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
A class for handling words, derived from string.
Definition: word.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual bool read(const dictionary &dict)
Read dictionary.
Namespace for OpenFOAM.
variableHeatTransfer(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from dictionary.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
addToRunTimeSelectionTable(option, fixedTemperatureConstraint, dictionary)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:49
virtual ~variableHeatTransfer()
Destructor.
labelList fv(nPoints)
Base class for inter region heat exchange. The derived classes must provide the heat transfer coeffis...
static const word propertiesName
Default name of the turbulence properties dictionary.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...