interRegionHeatTransferModel.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-2018 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 
27 #include "basicThermo.H"
28 #include "fvmSup.H"
30 #include "fvcVolumeIntegrate.H"
31 #include "fvOptionList.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace fv
38 {
39  defineTypeNameAndDebug(interRegionHeatTransferModel, 0);
40 }
41 }
42 
43 
44 // * * * * * * * * * * * * Protected member functions * * * * * * * * * * * //
45 
47 {
48  if (!firstIter_)
49  {
50  return;
51  }
52 
53  const fvMesh& nbrMesh = mesh_.time().lookupObject<fvMesh>(nbrRegionName_);
54 
55  const optionList& fvOptions = nbrMesh.lookupObject<optionList>("fvOptions");
56 
57  bool nbrModelFound = false;
58 
59  forAll(fvOptions, i)
60  {
61  if (fvOptions[i].name() == nbrModelName_)
62  {
63  nbrModel_ = &const_cast<interRegionHeatTransferModel&>
64  (
65  refCast<const interRegionHeatTransferModel>(fvOptions[i])
66  );
67  nbrModelFound = true;
68  break;
69  }
70  }
71 
72  if (!nbrModelFound)
73  {
75  << "Neighbour model not found" << nbrModelName_
76  << " in region " << nbrMesh.name() << nl
77  << exit(FatalError);
78  }
79 
80  firstIter_ = false;
81 
82  // Set nbr model's nbr model to avoid construction order problems
83  nbrModel_->setNbrModel();
84 }
85 
86 
88 {
89  if (master_)
90  {
91  if (mesh_.time().timeIndex() != timeIndex_)
92  {
93  calculateHtc();
94  timeIndex_ = mesh_.time().timeIndex();
95  }
96  }
97  else
98  {
99  nbrModel().correct();
100  interpolate(nbrModel().htc(), htc_);
101  }
102 }
103 
104 
105 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
106 
108 (
109  const word& name,
110  const word& modelType,
111  const dictionary& dict,
112  const fvMesh& mesh
113 )
114 :
116  (
117  name,
118  modelType,
119  dict,
120  mesh
121  ),
122  nbrModelName_(coeffs_.lookup("nbrModel")),
123  nbrModel_(nullptr),
124  firstIter_(true),
125  timeIndex_(-1),
126  htc_
127  (
128  IOobject
129  (
130  type() + ":htc",
131  mesh.time().timeName(),
132  mesh,
135  ),
136  mesh,
138  (
139  "htc",
141  0.0
142  ),
143  zeroGradientFvPatchScalarField::typeName
144  ),
145  semiImplicit_(false),
146  TName_(coeffs_.lookupOrDefault<word>("T", "T")),
147  TNbrName_(coeffs_.lookupOrDefault<word>("TNbr", "T"))
148 {
149  if (active())
150  {
151  coeffs_.lookup("fields") >> fieldNames_;
152  applied_.setSize(fieldNames_.size(), false);
153 
154  coeffs_.lookup("semiImplicit") >> semiImplicit_;
155  }
156 }
157 
158 
159 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
160 
162 {}
163 
164 
165 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
166 
168 (
169  fvMatrix<scalar>& eqn,
170  const label fieldi
171 )
172 {
173  setNbrModel();
174 
175  correct();
176 
177  const volScalarField& he = eqn.psi();
178 
179  const volScalarField& T = mesh_.lookupObject<volScalarField>(TName_);
180 
181  tmp<volScalarField> tTmapped
182  (
183  new volScalarField
184  (
185  IOobject
186  (
187  type() + ":Tmapped",
188  mesh_.time().timeName(),
189  mesh_,
192  ),
193  T
194  )
195  );
196 
197  volScalarField& Tmapped = tTmapped.ref();
198 
199  const fvMesh& nbrMesh = mesh_.time().lookupObject<fvMesh>(nbrRegionName_);
200 
201  const volScalarField& Tnbr =
202  nbrMesh.lookupObject<volScalarField>(TNbrName_);
203 
204  interpolate(Tnbr, Tmapped.primitiveFieldRef());
205 
206  if (debug)
207  {
208  Info<< "Volumetric integral of htc: "
209  << fvc::domainIntegrate(htc_).value()
210  << endl;
211 
212  if (mesh_.time().writeTime())
213  {
214  Tmapped.write();
215  htc_.write();
216  }
217  }
218 
219  if (semiImplicit_)
220  {
221  if (he.dimensions() == dimEnergy/dimMass)
222  {
223  if (mesh_.foundObject<basicThermo>(basicThermo::dictName))
224  {
225  const basicThermo& thermo =
226  mesh_.lookupObject<basicThermo>(basicThermo::dictName);
227 
228  volScalarField htcByCpv(htc_/thermo.Cpv());
229 
230  eqn += htc_*(Tmapped - T) + htcByCpv*he - fvm::Sp(htcByCpv, he);
231 
232  if (debug)
233  {
234  const dimensionedScalar energy =
235  fvc::domainIntegrate(htc_*(he/thermo.Cp() - Tmapped));
236 
237  Info<< "Energy exchange from region " << nbrMesh.name()
238  << " To " << mesh_.name() << " : " << energy.value()
239  << endl;
240  }
241  }
242  else
243  {
245  << " on mesh " << mesh_.name()
246  << " could not find object basicThermo."
247  << " The available objects are: "
248  << mesh_.names()
249  << exit(FatalError);
250  }
251  }
252  else if (he.dimensions() == dimTemperature)
253  {
254  eqn += htc_*Tmapped - fvm::Sp(htc_, he);
255  }
256  }
257  else
258  {
259  eqn += htc_*(Tmapped - T);
260  }
261 }
262 
263 
265 (
266  const volScalarField& rho,
267  fvMatrix<scalar>& eqn,
268  const label fieldi
269 )
270 {
271  addSup(eqn, fieldi);
272 }
273 
274 
275 // ************************************************************************* //
defineTypeNameAndDebug(option, 0)
volScalarField & he
Definition: YEEqn.H:51
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
Base class for inter region heat exchange. The derived classes must provide the heat transfer coeffis...
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
void setNbrModel()
Set the neighbour interRegionHeatTransferModel.
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
fv::options & fvOptions
Base class for inter-region exchange.
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:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
interRegionHeatTransferModel(const word &name, const word &modelType, const dictionary &dict, const fvMesh &mesh)
Construct from dictionary.
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
const GeometricField< Type, fvPatchField, volMesh > & psi() const
Definition: fvMatrix.H:284
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
dimensioned< Type > domainIntegrate(const GeometricField< Type, fvPatchField, volMesh > &vf)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
rhoReactionThermo & thermo
Definition: createFields.H:28
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
virtual tmp< volScalarField > Cpv() const =0
Heat capacity at constant pressure/volume [J/kg/K].
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Source term to energy equation.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:115
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
const dimensionSet & dimensions() const
Return dimensions.
dynamicFvMesh & mesh
void correct()
Correct to calculate the inter-region heat transfer coefficient.
A class for handling words, derived from string.
Definition: word.H:59
labelList fv(nPoints)
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
const Type & value() const
Return const reference to value.
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 const char nl
Definition: Ostream.H:265
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-6/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:72
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
virtual tmp< volScalarField > Cp() const =0
Heat capacity at constant pressure [J/kg/K].
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const dimensionSet dimEnergy
Internal & ref()
Return a reference to the dimensioned internal field.
const word & name() const
Return reference to name.
Definition: fvMesh.H:257
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
messageStream Info
virtual Ostream & write(const token &)=0
Write next token to stream.
List of finite volume options.
Definition: fvOptionList.H:64
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Calculate the matrix for implicit and explicit sources.
Namespace for OpenFOAM.