All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
interfaceCompositionModel.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-2020 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 "phaseModel.H"
28 #include "phasePair.H"
29 #include "phaseSystem.H"
30 #include "rhoReactionThermo.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(interfaceCompositionModel, 0);
37  defineRunTimeSelectionTable(interfaceCompositionModel, dictionary);
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const dictionary& dict,
46  const phasePair& pair
47 )
48 :
49  pair_(pair),
50  species_(dict.lookup("species")),
51  Le_("Le", dimless, dict),
52  phase_
53  (
54  // !!! This is a hack, to let the models know which side of the
55  // interface it applies to. In general, we are going to need improve
56  // the phase-pair system to take into account model sidedness and other
57  // types of asymmetry from just dispersed/continuous.
58  pair.phase1().name() == IOobject::group(dict.dictName())
59  ? pair.phase1()
60  : pair.phase2()
61  ),
62  otherPhase_(pair.otherPhase(phase_)),
63  thermo_
64  (
65  phase_.mesh().lookupObject<rhoReactionThermo>
66  (
67  IOobject::groupName(basicThermo::dictName, phase_.name())
68  )
69  ),
70  otherThermo_
71  (
72  otherPhase_.mesh().lookupObject<rhoThermo>
73  (
74  IOobject::groupName(basicThermo::dictName, otherPhase_.name())
75  )
76  )
77 {}
78 
79 
80 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
81 
83 {}
84 
85 
86 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
87 
89 (
90  const word& speciesName,
91  const volScalarField& Tf
92 ) const
93 {
94  const label speciei = composition().species()[speciesName];
95 
96  return Yf(speciesName, Tf) - composition().Y()[speciei];
97 }
98 
99 
101 (
102  const word& speciesName,
103  const volScalarField& Tf
104 ) const
105 {
106  return YfPrime(speciesName, Tf);
107 }
108 
109 
111 (
112  const word& speciesName
113 ) const
114 {
115  const label speciei = composition().species()[speciesName];
116  const volScalarField& p(thermo_.p());
117  const volScalarField& T(thermo_.T());
118 
119  return volScalarField::New
120  (
121  IOobject::groupName("D", pair_.name()),
122  composition().kappa(speciei, p, T)
123  /composition().Cp(speciei, p, T)
124  /composition().rho(speciei, p, T)
125  /Le_
126  );
127 }
128 
129 
130 // ************************************************************************* //
const char *const group
Group name for atomic constants.
dictionary dict
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
tmp< volScalarField > D(const word &speciesName) const
Mass diffusivity.
virtual scalar kappa(const label speciei, const scalar p, const scalar T) const =0
Thermal conductivity [W/m/K].
interfaceCompositionModel(const dictionary &dict, const phasePair &pair)
Construct from a dictionary and a phase pair.
tmp< volScalarField > dYfPrime(const word &speciesName, const volScalarField &Tf) const
Mass fraction difference between the interface and the field.
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
const dimensionSet dimless
virtual scalar rho(const label speciei, const scalar p, const scalar T) const =0
Density [kg/m^3].
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
virtual scalar Cp(const label speciei, const scalar p, const scalar T) const =0
Heat capacity at constant pressure [J/kg/K].
static word groupName(Name name, const word &group)
virtual ~interfaceCompositionModel()
Destructor.
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
defineTypeNameAndDebug(combustionModel, 0)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const basicSpecieMixture & composition() const
Return the composition.
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
tmp< volScalarField > dY(const word &speciesName, const volScalarField &Tf) const
Mass fraction difference between the interface and the field.
const word dictName("noiseDict")
Namespace for OpenFOAM.
const speciesTable & species() const
Return the table of species.