buoyancyEnergy.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) 2015-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 
26 #include "buoyancyEnergy.H"
27 #include "fvMatrices.H"
28 #include "basicThermo.H"
30 
31 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace fv
36 {
38 
40  (
41  fvModel,
44  );
45 }
46 }
47 
48 
49 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
50 
51 void Foam::fv::buoyancyEnergy::readCoeffs()
52 {
53  phaseName_ = coeffs().lookupOrDefault<word>("phase", word::null);
54 
55  UName_ =
56  coeffs().lookupOrDefault<word>
57  (
58  "U",
59  IOobject::groupName("U", phaseName_)
60  );
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
67 (
68  const word& name,
69  const word& modelType,
70  const fvMesh& mesh,
71  const dictionary& dict
72 )
73 :
74  fvModel(name, modelType, mesh, dict),
75  phaseName_(word::null),
76  UName_(word::null)
77 {
78  readCoeffs();
79 }
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
85 {
86  const basicThermo& thermo =
87  mesh().lookupObject<basicThermo>
88  (
89  IOobject::groupName(physicalProperties::typeName, phaseName_)
90  );
91 
92  return wordList(1, thermo.he().name());
93 }
94 
95 
97 (
98  const volScalarField& rho,
99  const volScalarField& he,
100  fvMatrix<scalar>& eqn
101 ) const
102 {
104  mesh().lookupObject<uniformDimensionedVectorField>("g");
105 
106  const volVectorField& U = mesh().lookupObject<volVectorField>(UName_);
107 
108  eqn += rho*(U&g);
109 }
110 
111 
113 (
114  const volScalarField& alpha,
115  const volScalarField& rho,
116  const volScalarField& he,
117  fvMatrix<scalar>& eqn
118 ) const
119 {
121  mesh().lookupObject<uniformDimensionedVectorField>("g");
122 
123  const volVectorField& U = mesh().lookupObject<volVectorField>(UName_);
124 
125  eqn += alpha*rho*(U&g);
126 }
127 
128 
130 {
131  return true;
132 }
133 
134 
136 {}
137 
138 
140 {}
141 
142 
144 {}
145 
146 
148 {
149  if (fvModel::read(dict))
150  {
151  readCoeffs();
152  return true;
153  }
154  else
155  {
156  return false;
157  }
158 }
159 
160 
161 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
static word groupName(Name name, const word &group)
Base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:78
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.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:118
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
Calculates and applies the buoyancy energy source rho*(U&g) to the energy equation.
virtual bool movePoints()
Update for mesh motion.
virtual wordList addSupFields() const
Return the list of fields for which the fvModel adds source term.
buoyancyEnergy(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict)
Construct from explicit source name and mesh.
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 source dictionary.
virtual void addSup(const volScalarField &rho, const volScalarField &he, fvMatrix< scalar > &eqn) const
Add explicit contribution to compressible energy equation.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
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.
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
A special matrix type and solver, designed for finite volume solutions of scalar equations.
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
addToRunTimeSelectionTable(fvConstraint, bound, dictionary)
defineTypeNameAndDebug(bound, 0)
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
thermo he()
labelList fv(nPoints)
dictionary dict
fluidMulticomponentThermo & thermo
Definition: createFields.H:31