heatTransferModel.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) 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 
26 #include "heatTransferModel.H"
27 
28 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace fv
33 {
34  defineTypeNameAndDebug(heatTransferModel, 0);
35  defineRunTimeSelectionTable(heatTransferModel, mesh);
36  defineRunTimeSelectionTable(heatTransferModel, model);
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 void Foam::fv::heatTransferModel::readCoeffs()
44 {
45  IOobject AoVIO
46  (
47  "AoV",
48  mesh().time().constant(),
49  mesh(),
52  );
53 
54  if (coeffs().found("AoV"))
55  {
56  AoV_ = dimensionedScalar("AoV", dimless/dimLength, coeffs());
57  AoVPtr_.clear();
58  }
59  else if (AoVIO.typeHeaderOk<volScalarField>(false))
60  {
61  AoV_ = dimensionedScalar("AoV", dimless/dimLength, NaN);
62  AoVPtr_.set(new volScalarField(AoVIO, mesh()));
63  }
64  else
65  {
66  FatalIOErrorInFunction(coeffs())
67  << "Area per unit volume (AoV) not found. A uniform AoV "
68  << "value should be specified, or a non-uniform field should "
69  << "exist at " << AoVIO.objectPath()
70  << exit(FatalIOError);
71  }
72 }
73 
74 
75 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
76 
78 (
79  const word& modelType,
80  const dictionary& dict,
81  const fvMesh& mesh
82 )
83 :
84  mesh_(mesh),
85  coeffs_(dict.optionalSubDict(modelType + "Coeffs")),
86  AoV_("AoV", dimless/dimLength, NaN),
87  AoVPtr_(nullptr)
88 {
89  readCoeffs();
90 }
91 
92 
94 (
95  const word& modelType,
96  const dictionary& dict,
97  const interRegionModel& model
98 )
99 :
100  heatTransferModel(modelType, dict, model.mesh())
101 {}
102 
103 
104 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
105 
108 {
109  word heatTransferModelType(dict.lookup("type"));
110 
111  Info<< "Selecting heatTransferModel "
112  << heatTransferModelType << endl;
113 
114  meshConstructorTable::iterator cstrIter =
115  meshConstructorTablePtr_->find(heatTransferModelType);
116 
117  if (cstrIter == meshConstructorTablePtr_->end())
118  {
120  << "Unknown heatTransferModelType type "
121  << heatTransferModelType << endl << endl
122  << "Valid heatTransferModel types are : " << endl
123  << meshConstructorTablePtr_->sortedToc()
124  << exit(FatalError);
125  }
126 
127  return cstrIter()(dict, mesh);
128 }
129 
130 
133 (
134  const dictionary& dict,
135  const interRegionModel& model
136 )
137 {
138  word heatTransferModelType(dict.lookup("type"));
139 
140  Info<< "Selecting heatTransferModel "
141  << heatTransferModelType << endl;
142 
143  modelConstructorTable::iterator cstrIter =
144  modelConstructorTablePtr_->find(heatTransferModelType);
145 
146  if (cstrIter == modelConstructorTablePtr_->end())
147  {
149  << "Unknown heatTransferModelType type "
150  << heatTransferModelType << endl << endl
151  << "Valid heatTransferModel types are : " << endl
152  << modelConstructorTablePtr_->sortedToc()
153  << exit(FatalError);
154  }
155 
156  return cstrIter()(dict, model);
157 }
158 
159 
160 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
161 
163 {}
164 
165 
166 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
167 
169 {
170  if (!AoVPtr_.valid())
171  {
172  return volScalarField::New(type() + ":AoV", mesh_, AoV_);
173  }
174  else
175  {
176  return AoVPtr_();
177  }
178 }
179 
180 
182 {
183  coeffs_ = dict.optionalSubDict(type() + "Coeffs");
184 
185  readCoeffs();
186 
187  return true;
188 }
189 
190 
191 // ************************************************************************* //
tmp< volScalarField > AoV() const
Get the area per unit volume.
dictionary dict
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
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
Base class for heat transfer coefficient modelling used in heat transfer fvModels. Area per unit volume [1/m] (AoV) must be provided as a value in the coefficients dictionary or as a field in constant.
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
defineRunTimeSelectionTable(heatTransferModel, mesh)
const dimensionSet dimless
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvModelI.H:34
const dimensionSet dimLength
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Base class for inter-region exchange.
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1040
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
static autoPtr< heatTransferModel > New(const dictionary &dict, const fvMesh &mesh)
Select from dictionary and mesh.
heatTransferModel(const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from dictionary and mesh.
virtual bool read(const dictionary &dict)
Read dictionary.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:335
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
messageStream Info
A class for managing temporary objects.
Definition: PtrList.H:53
bool found
virtual ~heatTransferModel()
Destructor.
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
IOerror FatalIOError