heSolidThermo.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 
26 #include "heSolidThermo.H"
27 #include "volFields.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class BasicSolidThermo, class MixtureType>
33 {
34  scalarField& TCells = this->T_.primitiveFieldRef();
35 
36  const scalarField& hCells = this->he_;
37  const scalarField& pCells = this->p_;
38  scalarField& rhoCells = this->rho_.primitiveFieldRef();
39  scalarField& alphaCells = this->alpha_.primitiveFieldRef();
40 
41  forAll(TCells, celli)
42  {
43  const typename MixtureType::thermoType& mixture_ =
44  this->cellMixture(celli);
45 
46  const typename MixtureType::thermoType& volMixture_ =
47  this->cellVolMixture(pCells[celli], TCells[celli], celli);
48 
49  TCells[celli] = mixture_.THE
50  (
51  hCells[celli],
52  pCells[celli],
53  TCells[celli]
54  );
55 
56  rhoCells[celli] = volMixture_.rho(pCells[celli], TCells[celli]);
57 
58  alphaCells[celli] =
59  volMixture_.kappa(pCells[celli], TCells[celli])
60  /
61  mixture_.Cpv(pCells[celli], TCells[celli]);
62  }
63 
64  volScalarField::Boundary& pBf =
65  this->p_.boundaryFieldRef();
66 
67  volScalarField::Boundary& TBf =
68  this->T_.boundaryFieldRef();
69 
70  volScalarField::Boundary& rhoBf =
71  this->rho_.boundaryFieldRef();
72 
73  volScalarField::Boundary& heBf =
74  this->he().boundaryFieldRef();
75 
76  volScalarField::Boundary& alphaBf =
77  this->alpha_.boundaryFieldRef();
78 
79  forAll(this->T_.boundaryField(), patchi)
80  {
81  fvPatchScalarField& pp = pBf[patchi];
82  fvPatchScalarField& pT = TBf[patchi];
83  fvPatchScalarField& prho = rhoBf[patchi];
84  fvPatchScalarField& phe = heBf[patchi];
85  fvPatchScalarField& palpha = alphaBf[patchi];
86 
87  if (pT.fixesValue())
88  {
89  forAll(pT, facei)
90  {
91  const typename MixtureType::thermoType& mixture_ =
92  this->patchFaceMixture(patchi, facei);
93 
94  const typename MixtureType::thermoType& volMixture_ =
95  this->patchFaceVolMixture
96  (
97  pp[facei],
98  pT[facei],
99  patchi,
100  facei
101  );
102 
103 
104  phe[facei] = mixture_.HE(pp[facei], pT[facei]);
105  prho[facei] = volMixture_.rho(pp[facei], pT[facei]);
106 
107  palpha[facei] =
108  volMixture_.kappa(pp[facei], pT[facei])
109  / mixture_.Cpv(pp[facei], pT[facei]);
110  }
111  }
112  else
113  {
114  forAll(pT, facei)
115  {
116  const typename MixtureType::thermoType& mixture_ =
117  this->patchFaceMixture(patchi, facei);
118 
119  const typename MixtureType::thermoType& volMixture_ =
120  this->patchFaceVolMixture
121  (
122  pp[facei],
123  pT[facei],
124  patchi,
125  facei
126  );
127 
128  pT[facei] = mixture_.THE(phe[facei], pp[facei] ,pT[facei]);
129  prho[facei] = volMixture_.rho(pp[facei], pT[facei]);
130 
131  palpha[facei] =
132  volMixture_.kappa(pp[facei], pT[facei])
133  / mixture_.Cpv(pp[facei], pT[facei]);
134  }
135  }
136  }
137 
138  this->alpha_.correctBoundaryConditions();
139 }
140 
141 
142 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
143 
144 template<class BasicSolidThermo, class MixtureType>
147 (
148  const fvMesh& mesh,
149  const word& phaseName
150 )
151 :
153 {
154  calculate();
155 }
156 
157 
158 template<class BasicSolidThermo, class MixtureType>
161 (
162  const fvMesh& mesh,
163  const dictionary& dict,
164  const word& phaseName
165 )
166 :
168 {
169  calculate();
170 }
171 
172 
173 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
174 
175 template<class BasicSolidThermo, class MixtureType>
177 {}
178 
179 
180 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
181 
182 template<class BasicSolidThermo, class MixtureType>
184 {
185  if (debug)
186  {
187  InfoInFunction << endl;
188  }
189 
190  calculate();
191 
192  if (debug)
193  {
194  Info<< " Finished" << endl;
195  }
196 }
197 
198 
199 template<class BasicSolidThermo, class MixtureType>
202 {
203  const fvMesh& mesh = this->T_.mesh();
204 
205  tmp<volVectorField> tKappa
206  (
207  new volVectorField
208  (
209  IOobject
210  (
211  "Kappa",
212  mesh.time().timeName(),
213  mesh,
214  IOobject::NO_READ,
215  IOobject::NO_WRITE
216  ),
217  mesh,
219  )
220  );
221 
222  volVectorField& Kappa = tKappa.ref();
223  vectorField& KappaCells = Kappa.primitiveFieldRef();
224  const scalarField& TCells = this->T_;
225  const scalarField& pCells = this->p_;
226 
227  forAll(KappaCells, celli)
228  {
229  Kappa[celli] =
230  this->cellVolMixture
231  (
232  pCells[celli],
233  TCells[celli],
234  celli
235  ).Kappa(pCells[celli], TCells[celli]);
236  }
237 
238  volVectorField::Boundary& KappaBf = Kappa.boundaryFieldRef();
239 
240  forAll(KappaBf, patchi)
241  {
242  vectorField& Kappap = KappaBf[patchi];
243  const scalarField& pT = this->T_.boundaryField()[patchi];
244  const scalarField& pp = this->p_.boundaryField()[patchi];
245 
246  forAll(Kappap, facei)
247  {
248  Kappap[facei] =
249  this->patchFaceVolMixture
250  (
251  pp[facei],
252  pT[facei],
253  patchi,
254  facei
255  ).Kappa(pp[facei], pT[facei]);
256  }
257  }
258 
259  return tKappa;
260 }
261 
262 
263 template<class BasicSolidThermo, class MixtureType>
266 (
267  const label patchi
268 ) const
269 {
270  const scalarField& pp = this->p_.boundaryField()[patchi];
271  const scalarField& Tp = this->T_.boundaryField()[patchi];
272  tmp<vectorField> tKappa(new vectorField(pp.size()));
273 
274  vectorField& Kappap = tKappa.ref();
275 
276  forAll(Tp, facei)
277  {
278  Kappap[facei] =
279  this->patchFaceVolMixture
280  (
281  pp[facei],
282  Tp[facei],
283  patchi,
284  facei
285  ).Kappa(pp[facei], Tp[facei]);
286  }
287 
288  return tKappa;
289 }
290 
291 
292 // ************************************************************************* //
volScalarField & he
Definition: YEEqn.H:51
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
volVectorField vectorField(fieldObject, mesh)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
virtual void correct()
Update properties.
virtual tmp< volVectorField > Kappa() const
Anisotropic thermal conductivity [W/m/K].
dynamicFvMesh & mesh
fvPatchField< scalar > fvPatchScalarField
A class for handling words, derived from string.
Definition: word.H:59
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:52
volScalarField scalarField(fieldObject, mesh)
const dimensionSet dimEnergy
Internal & ref()
Return a reference to the dimensioned internal field.
label patchi
Energy for a solid mixture.
Definition: heSolidThermo.H:49
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
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
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
virtual ~heSolidThermo()
Destructor.
#define InfoInFunction
Report an information message using Foam::Info.