singleComponentPhaseChange.H
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-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 Class
25  Foam::fv::singleComponentPhaseChange
26 
27 Description
28  Base class for phase change models in which only a single component changes
29  phase. Can be applied to any combination of pure and multicomponent phases.
30  If either phase is multicomponent, then a single specie must be identified
31  as the one that changes phase.
32 
33 SourceFiles
34  singleComponentPhaseChange.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef singleComponentPhaseChange_H
39 #define singleComponentPhaseChange_H
40 
41 #include "phaseChange.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 namespace fv
48 {
49 
50 /*---------------------------------------------------------------------------*\
51  Class singleComponentPhaseChange Declaration
52 \*---------------------------------------------------------------------------*/
53 
55 :
56  public phaseChange
57 {
58 private:
59 
60  // Private Data
61 
62  //- The name of the changing specie, or word::null if neither thermo
63  // is multicomponent
64  const word specie_;
65 
66  //- The indices of the changing species in the phases, or -1 if the
67  // phases are not multicomponent
68  const Pair<label> specieis_;
69 
70  //- Whether or not to linearise the energy source
71  bool energySemiImplicit_;
72 
73 
74  // Private Member Functions
75 
76  //- Non-virtual read
77  void readCoeffs();
78 
79 
80 public:
81 
82  //- Runtime type information
83  TypeName("singleComponentPhaseChange");
84 
85 
86  // Constructors
87 
88  //- Construct from explicit source name and mesh
90  (
91  const word& name,
92  const word& modelType,
93  const fvMesh& mesh,
94  const dictionary& dict,
95  const Pair<bool>& fluidThermosRequired,
96  const Pair<bool>& specieThermosRequired
97  );
98 
99 
100  // Member Functions
101 
102  // Access
103 
104  //- Return the name of the changing specie
105  inline const word& specie() const;
106 
107  //- Return the indices of the changing species in the phases
108  inline const Pair<label>& specieis() const;
109 
110 
111  // Sources
112 
113  //- Override the energy equation to add the phase change heat, or
114  // the species equations to add the relevant mass sources
115  void addSup
116  (
117  const volScalarField& alpha,
118  const volScalarField& rho,
119  const volScalarField& heOrYi,
120  fvMatrix<scalar>& eqn
121  ) const;
122 
123 
124  //- Read source dictionary
125  virtual bool read(const dictionary& dict);
126 };
127 
128 
129 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
130 
131 } // End namespace fv
132 } // End namespace Foam
133 
134 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
135 
137 
138 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
139 
140 #endif
141 
142 // ************************************************************************* //
Generic GeometricField class.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
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
const fvMesh & mesh() const
Return const access to the mesh database.
Definition: fvModelI.H:53
const word & name() const
Return const access to the source name.
Definition: fvModelI.H:47
tmp< volScalarField::Internal > rho(const label i) const
Return the density.
Definition: massTransfer.C:85
Base class for phase change models.
Definition: phaseChange.H:59
Base class for phase change models in which only a single component changes phase....
const Pair< label > & specieis() const
Return the indices of the changing species in the phases.
const word & specie() const
Return the name of the changing specie.
virtual bool read(const dictionary &dict)
Read source dictionary.
void addSup(const volScalarField &alpha, const volScalarField &rho, const volScalarField &heOrYi, fvMatrix< scalar > &eqn) const
Override the energy equation to add the phase change heat, or.
TypeName("singleComponentPhaseChange")
Runtime type information.
singleComponentPhaseChange(const word &name, const word &modelType, const fvMesh &mesh, const dictionary &dict, const Pair< bool > &fluidThermosRequired, const Pair< bool > &specieThermosRequired)
Construct from explicit source name and mesh.
A class for handling words, derived from string.
Definition: word.H:62
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Namespace for OpenFOAM.
labelList fv(nPoints)
dictionary dict