IsothermalPhaseModel.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-2023 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 "IsothermalPhaseModel.H"
27 #include "phaseSystem.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class BasePhaseModel>
33 (
34  const phaseSystem& fluid,
35  const word& phaseName,
36  const bool referencePhase,
37  const label index
38 )
39 :
40  BasePhaseModel(fluid, phaseName, referencePhase, index)
41 {}
42 
43 
44 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
45 
46 template<class BasePhaseModel>
48 {}
49 
50 
51 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
52 
53 template<class BasePhaseModel>
55 {
56  BasePhaseModel::correctThermo();
57 
58  // Correct the thermo for pressure changes
59  // ensuring the temperature remains constant
61  (
63  (
64  this->thermo().T().name() + ":Copy",
65  this->thermo().T()
66  )
67  );
68  this->thermo_->he() = this->thermo().he(this->fluidThermo().p(), TCopy);
69  this->thermo_->correct();
70  this->thermo_->T() = TCopy;
71 }
72 
73 
74 template<class BasePhaseModel>
76 {
77  return true;
78 }
79 
80 
81 template<class BasePhaseModel>
84 {
86  return this->thermo().kappa().boundaryField()[patchi];
87 }
88 
89 
90 template<class BasePhaseModel>
93 {
95  << "Cannot construct an energy equation for an isothermal phase"
96  << exit(FatalError);
97 
98  return tmp<fvScalarMatrix>();
99 }
100 
101 
102 // ************************************************************************* //
tmp< GeometricField< Type, PatchField, GeoMesh > > T() const
Return transpose (only if it is a tensor field)
const Boundary & boundaryField() const
Return const-reference to the boundary field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
virtual void correctThermo()
Correct the thermodynamics.
virtual ~IsothermalPhaseModel()
Destructor.
IsothermalPhaseModel(const phaseSystem &fluid, const word &phaseName, const bool referencePhase, const label index)
virtual bool isothermal() const
Return whether the phase is isothermal.
virtual tmp< fvScalarMatrix > heEqn()
Return the enthalpy equation.
virtual tmp< scalarField > kappaEff(const label patchi) const
Return the effective thermal conductivity on a patch.
virtual const volScalarField & he() const =0
Enthalpy/Internal energy [J/kg].
virtual const volScalarField & kappa() const =0
Thermal conductivity of mixture [W/m/K].
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:57
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:73
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
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
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
error FatalError
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
volScalarField & p
fluidMulticomponentThermo & thermo
Definition: createFields.H:31