All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Henry.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 
26 #include "Henry.H"
27 #include "phasePair.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace interfaceCompositionModels
35 {
36  defineTypeNameAndDebug(Henry, 0);
37  addToRunTimeSelectionTable(interfaceCompositionModel, Henry, dictionary);
38 }
39 }
40 
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
45 (
46  const dictionary& dict,
47  const phasePair& pair
48 )
49 :
50  interfaceCompositionModel(dict, pair),
51  k_(dict.lookup("k")),
52  YSolvent_
53  (
54  IOobject
55  (
56  IOobject::groupName("YSolvent", pair.name()),
57  pair.phase1().mesh().time().timeName(),
58  pair.phase1().mesh()
59  ),
60  pair.phase1().mesh(),
62  )
63 {
64  if (k_.size() != species().size())
65  {
67  << "Differing number of species and solubilities"
68  << exit(FatalError);
69  }
70 }
71 
72 
73 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
74 
76 {}
77 
78 
79 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
80 
82 {
83  YSolvent_ = scalar(1);
84 
85  forAllConstIter(hashedWordList, species(), iter)
86  {
87  YSolvent_ -= Yf(*iter, Tf);
88  }
89 }
90 
91 
93 (
94  const word& speciesName,
95  const volScalarField& Tf
96 ) const
97 {
98  if (species().found(speciesName))
99  {
100  const label index = species()[speciesName];
101 
102  return
103  k_[index]
104  *otherComposition().Y(speciesName)
105  *otherThermo().rho()
106  /thermo().rho();
107  }
108  else
109  {
110  return YSolvent_*composition().Y(speciesName);
111  }
112 }
113 
114 
116 (
117  const word& speciesName,
118  const volScalarField& Tf
119 ) const
120 {
121  return volScalarField::New
122  (
123  IOobject::groupName("YfPrime", pair().name()),
124  pair().phase1().mesh(),
126  );
127 }
128 
129 
130 // ************************************************************************* //
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
fluidReactionThermo & thermo
Definition: createFields.H:28
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
basicSpecieMixture & composition
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 tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
static word groupName(Name name, const word &group)
Henry(const dictionary &dict, const phasePair &pair)
Construct from components.
word timeName
Definition: getTimeIndex.H:3
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void update(const volScalarField &Tf)
Update the composition.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
A class for managing temporary objects.
Definition: PtrList.H:53
bool found
const dimensionSet dimTemperature
Namespace for OpenFOAM.