solidThermalEquilibrium.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) 2019-2024 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 "fvmDdt.H"
28 #include "fvmLaplacian.H"
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace fv
36 {
40  (
41  fvModel,
43  dictionary,
44  solidEquilibriumEnergySource,
45  "solidEquilibriumEnergySource"
46  );
47 }
48 }
49 
50 
51 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
52 
53 void Foam::fv::solidThermalEquilibrium::readCoeffs()
54 {
55  phaseName_ = coeffs().lookupOrDefault<word>("phase", word::null);
56 
57  solidPhaseName_ = coeffs().lookup<word>("solidPhase");
58 }
59 
60 
62 Foam::fv::solidThermalEquilibrium::solidAlpha() const
63 {
64  const word alphaName = IOobject::groupName("alpha", solidPhaseName_);
65 
66  if (!mesh().foundObject<volScalarField>(alphaName))
67  {
68  volScalarField* alphaPtr =
69  new volScalarField
70  (
71  IOobject
72  (
73  alphaName,
74  mesh().time().constant(),
75  mesh(),
78  ),
79  mesh()
80  );
81 
82  alphaPtr->store();
83  }
84 
85  return mesh().lookupObject<volScalarField>(alphaName);
86 }
87 
88 
89 const Foam::solidThermo&
90 Foam::fv::solidThermalEquilibrium::solidThermo() const
91 {
92  const word thermoName =
93  IOobject::groupName(physicalProperties::typeName, solidPhaseName_);
94 
95  if (!mesh().foundObject<Foam::solidThermo>(thermoName))
96  {
97  Foam::solidThermo* thermoPtr =
98  solidThermo::New(mesh(), solidPhaseName_).ptr();
99 
100  thermoPtr->properties().store();
101  }
102 
103  return mesh().lookupObject<Foam::solidThermo>(thermoName);
104 }
105 
106 
107 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
108 
110 (
111  const word& name,
112  const word& modelType,
113  const fvMesh& mesh,
114  const dictionary& dict
115 )
116 :
117  fvModel(name, modelType, mesh, dict),
118  phaseName_(word::null),
119  solidPhaseName_(word::null)
120 {
121  read(dict);
122  solidAlpha();
123  solidThermo();
124 }
125 
126 
127 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
128 
130 {}
131 
132 
133 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
134 
136 {
137  const basicThermo& thermo =
138  mesh().lookupObject<basicThermo>
139  (
140  IOobject::groupName(physicalProperties::typeName, phaseName_)
141  );
142 
143  return wordList(1, thermo.he().name());
144 }
145 
146 
148 (
149  const volScalarField& rho,
150  const volScalarField& he,
151  fvMatrix<scalar>& eqn
152 ) const
153 {
154  const volScalarField alphahe
155  (
156  "alphahe",
157  solidThermo().kappa()/solidThermo().Cpv()
158  );
159 
160  const volScalarField& A = solidAlpha();
161  const volScalarField B(1 - A);
162 
163  eqn -=
164  A/B*fvm::ddt(solidThermo().rho(), eqn.psi());
165  - 1/B*fvm::laplacian
166  (
167  A*alphahe,
168  eqn.psi(),
169  "laplacian(" + alphahe.name() + "," + eqn.psi().name() + ")"
170  );
171 }
172 
173 
175 (
176  const volScalarField& alpha,
177  const volScalarField& rho,
178  const volScalarField& he,
179  fvMatrix<scalar>& eqn
180 ) const
181 {
182  const volScalarField alphahe
183  (
184  "alphahe",
185  solidThermo().kappa()/solidThermo().Cpv()
186  );
187 
188  const volScalarField& A = solidAlpha();
189  const volScalarField B(1 - A);
190 
191  eqn -=
192  A/B*fvm::ddt(alpha, solidThermo().rho(), eqn.psi());
193  - 1/B*fvm::laplacian
194  (
195  A*alphahe,
196  eqn.psi(),
197  "laplacian(" + alphahe.name() + "," + eqn.psi().name() + ")"
198  );
199 }
200 
201 
203 {
204  return true;
205 }
206 
207 
209 (
210  const polyTopoChangeMap&
211 )
212 {}
213 
214 
216 {}
217 
218 
220 (
221  const polyDistributionMap&
222 )
223 {}
224 
225 
227 {
228  if (fvModel::read(dict))
229  {
230  readCoeffs();
231  return true;
232  }
233  else
234  {
235  return false;
236  }
237 }
238 
239 
240 // ************************************************************************* //
static const Foam::dimensionedScalar A("A", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("B", Foam::dimless, 18.678)
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
const word & name() const
Return name.
Definition: IOobject.H:310
static word groupName(Name name, const word &group)
Base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:78
virtual const IOdictionary & properties() const =0
Properties dictionary.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:710
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
VolField< Type > & psi()
Definition: fvMatrix.H:289
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
Finite volume model abstract base class.
Definition: fvModel.H:59
const dictionary & coeffs() const
Return dictionary.
Definition: fvModelI.H:59
virtual bool read(const dictionary &dict)
Read source dictionary.
Definition: fvModel.C:199
This fvModel adds the thermal inertia of a solid phase into the energy equation. It assumes that the ...
virtual bool movePoints()
Update for mesh motion.
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
virtual void topoChange(const polyTopoChangeMap &)
Update topology using the given map.
virtual void distribute(const polyDistributionMap &)
Redistribute or update using the given distribution map.
virtual bool read(const dictionary &dict)
Read dictionary.
virtual void addSup(const volScalarField &rho, const volScalarField &he, fvMatrix< scalar > &eqn) const
Explicit and implicit sources for compressible equations.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
solidThermalEquilibrium(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from components.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
void store()
Transfer ownership of this object to its registry.
Definition: regIOobjectI.H:40
Base-class for solid thermodynamic properties.
Definition: solidThermo.H:59
static autoPtr< solidThermo > New(const fvMesh &, const word &phaseName=word::null)
Standard selection based on fvMesh.
Definition: solidThermo.C:67
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
Calculate the matrix for the first temporal derivative.
Calculate the matrix for the laplacian of the field.
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
addToRunTimeSelectionTable(fvConstraint, bound, dictionary)
defineTypeNameAndDebug(bound, 0)
addBackwardCompatibleToRunTimeSelectionTable(fvConstraint, fixedTemperature, dictionary, fixedTemperatureConstraint, "fixedTemperatureConstraint")
tmp< fvMatrix< Type > > laplacian(const VolField< Type > &vf, const word &name)
Definition: fvmLaplacian.C:47
tmp< fvMatrix< Type > > ddt(const VolField< Type > &vf)
Definition: fvmDdt.C:46
Namespace for OpenFOAM.
List< word > wordList
A List of words.
Definition: fileName.H:54
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
thermo he()
labelList fv(nPoints)
dictionary dict
fluidMulticomponentThermo & thermo
Definition: createFields.H:31