heatTransfer.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 "heatTransfer.H"
27 #include "basicThermo.H"
28 #include "fvmSup.H"
30 #include "fvcVolumeIntegrate.H"
31 #include "fvModels.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace fv
39 {
40  defineTypeNameAndDebug(heatTransfer, 0);
42  (
43  fvModel,
44  heatTransfer,
45  dictionary
46  );
47 }
48 }
49 
50 
51 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52 
53 void Foam::fv::heatTransfer::readCoeffs()
54 {
55  semiImplicit_ = coeffs().lookup<bool>("semiImplicit");
56 
57  TName_ = coeffs().lookupOrDefault<word>("T", "T");
58 
59  Ta_ = dimensionedScalar("Ta", dimTemperature, coeffs());
60 
61  heatTransferModel_ = heatTransferModel::New(coeffs(), mesh());
62 }
63 
64 
65 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
66 
68 (
69  const word& name,
70  const word& modelType,
71  const dictionary& dict,
72  const fvMesh& mesh
73 )
74 :
75  fvModel(name, modelType, dict, mesh),
76  set_(coeffs(), mesh),
77  semiImplicit_(false),
78  TName_(word::null),
79  Ta_("Ta", dimTemperature, NaN),
80  heatTransferModel_(nullptr)
81 {
82  readCoeffs();
83 }
84 
85 
86 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
87 
89 {}
90 
91 
92 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
93 
95 {
96  const basicThermo& thermo =
97  mesh().lookupObject<basicThermo>(basicThermo::dictName);
98 
99  return wordList(1, thermo.he().name());
100 }
101 
102 
104 (
105  fvMatrix<scalar>& eqn,
106  const word& fieldName
107 ) const
108 {
109  const volScalarField& he = eqn.psi();
110 
111  const volScalarField& T =
112  mesh().lookupObject<volScalarField>(TName_);
113 
114  tmp<volScalarField> mask =
116  UIndirectList<scalar>(mask.ref().primitiveFieldRef(), set_.cells()) = 1;
117  const volScalarField htcAoV
118  (
119  mask*heatTransferModel_->htc()*heatTransferModel_->AoV()
120  );
121 
122  if (semiImplicit_)
123  {
124  if (he.dimensions() == dimEnergy/dimMass)
125  {
126  const basicThermo& thermo =
127  mesh().lookupObject<basicThermo>(basicThermo::dictName);
128 
129  const volScalarField htcAoVByCpv(htcAoV/thermo.Cpv());
130 
131  eqn += htcAoV*(Ta_ - T) + htcAoVByCpv*he - fvm::Sp(htcAoVByCpv, he);
132  }
133  else if (he.dimensions() == dimTemperature)
134  {
135  eqn += htcAoV*Ta_ - fvm::Sp(htcAoV, he);
136  }
137  }
138  else
139  {
140  eqn += htcAoV*(Ta_ - T);
141  }
142 }
143 
144 
146 (
147  const volScalarField& rho,
148  fvMatrix<scalar>& eqn,
149  const word& fieldName
150 ) const
151 {
152  addSup(eqn, fieldName);
153 }
154 
155 
157 {
158  heatTransferModel_->correct();
159 }
160 
161 
163 {
164  if (fvModel::read(dict))
165  {
166  readCoeffs();
167  return true;
168  }
169  else
170  {
171  return false;
172  }
173 }
174 
175 
176 // ************************************************************************* //
heatTransfer(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from dictionary.
Definition: heatTransfer.C:68
defineTypeNameAndDebug(fixedTemperatureConstraint, 0)
const word & name() const
Return name.
Definition: IOobject.H:303
fluidReactionThermo & thermo
Definition: createFields.H:28
Base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:77
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:282
virtual bool read(const dictionary &dict)
Read dictionary.
Definition: heatTransfer.C:162
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:158
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Finite volume model abstract base class.
Definition: fvModel.H:55
const dimensionSet dimless
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
virtual tmp< volScalarField > Cpv() const =0
Heat capacity at constant pressure/volume [J/kg/K].
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
virtual ~heatTransfer()
Destructor.
Definition: heatTransfer.C:88
static const word null
An empty word.
Definition: word.H:77
Volume integrate volField creating a volField.
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
static autoPtr< heatTransferModel > New(const dictionary &dict, const fvMesh &mesh)
Select from dictionary and mesh.
static const word dictName
Name of the thermophysical properties dictionary.
Definition: basicThermo.H:129
virtual void correct()
Correct the model.
Definition: heatTransfer.C:156
thermo he()
const dimensionSet dimEnergy
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const dimensionSet dimMass
addToRunTimeSelectionTable(fvConstraint, fixedTemperatureConstraint, dictionary)
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
Definition: heatTransfer.C:94
List< word > wordList
A List of words.
Definition: fileName.H:54
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void addSup(fvMatrix< scalar > &eqn, const word &fieldName) const
Source term to energy equation.
Definition: heatTransfer.C:104
Calculate the matrix for implicit and explicit sources.
const dimensionSet dimTemperature
Namespace for OpenFOAM.