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-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 "Henry.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace interfaceCompositionModels
34 {
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const dictionary& dict,
46  const phaseInterface& interface
47 )
48 :
49  interfaceCompositionModel(dict, interface),
50  k_(dict.lookup("k")),
51  YSolvent_
52  (
53  IOobject
54  (
55  IOobject::groupName("YSolvent", this->interface().name()),
56  interface.mesh().time().name(),
57  interface.mesh()
58  ),
59  interface.mesh(),
61  )
62 {
63  if (k_.size() != species().size())
64  {
66  << "Differing number of species and solubilities"
67  << exit(FatalError);
68  }
69 }
70 
71 
72 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
73 
75 {}
76 
77 
78 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
79 
81 {
82  YSolvent_ = scalar(1);
83 
84  forAllConstIter(hashedWordList, species(), iter)
85  {
86  YSolvent_ -= Yf(*iter, Tf);
87  }
88 }
89 
90 
92 (
93  const word& speciesName,
94  const volScalarField& Tf
95 ) const
96 {
97  if (species().found(speciesName))
98  {
99  const label index = species()[speciesName];
100 
101  return
102  k_[index]
103  *otherMulticomponentThermo().Y(speciesName)
104  *otherThermo().rho()/thermo().rho();
105  }
106  else
107  {
108  return YSolvent_*thermo().Y(speciesName);
109  }
110 }
111 
112 
114 (
115  const word& speciesName,
116  const volScalarField& Tf
117 ) const
118 {
119  return volScalarField::New
120  (
121  IOobject::groupName("YfPrime", interface().name()),
122  interface().mesh(),
124  );
125 }
126 
127 
128 // ************************************************************************* //
bool found
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
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,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static word groupName(Name name, const word &group)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
A wordList with hashed indices for faster lookup by name.
Generic base class for interface composition models. These models describe the composition in phase 1...
const hashedWordList & species() const
Return the transferring species names.
Henry's law for gas solubility in liquid. The concentration of a dissolved species in the liquid is p...
Definition: Henry.H:58
Henry(const dictionary &dict, const phaseInterface &interface)
Construct from a dictionary and an interface.
Definition: Henry.C:44
virtual ~Henry()
Destructor.
Definition: Henry.C:74
virtual void update(const volScalarField &Tf)
Update the composition.
Definition: Henry.C:80
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
Definition: Henry.C:92
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
Definition: Henry.C:114
Class to represent an interface between phases. Derivations can further specify the configuration of ...
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
addToRunTimeSelectionTable(interfaceCompositionModel, Henry, dictionary)
Namespace for OpenFOAM.
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
const dimensionSet dimless
const dimensionSet dimTemperature
error FatalError
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dictionary dict
fluidMulticomponentThermo & thermo
Definition: createFields.H:31