Henry.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015-2018 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "Henry.H"
29
30// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31
32template<class Thermo, class OtherThermo>
34(
35 const dictionary& dict,
36 const phasePair& pair
37)
38:
39 InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
40 k_(dict.lookup("k")),
41 YSolvent_
42 (
44 (
45 IOobject::groupName("YSolvent", pair.name()),
46 pair.phase1().mesh().time().timeName(),
47 pair.phase1().mesh()
48 ),
49 pair.phase1().mesh(),
50 dimensionedScalar("one", dimless, 1)
51 )
52{
53 if (k_.size() != this->speciesNames_.size())
54 {
56 << "Differing number of species and solubilities"
57 << exit(FatalError);
58 }
59}
60
61
62// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
63
64template<class Thermo, class OtherThermo>
66(
67 const volScalarField& Tf
68)
69{
70 YSolvent_ = scalar(1);
71
72 for (const word& speciesName : this->speciesNames_)
73 {
74 YSolvent_ -= Yf(speciesName, Tf);
75 }
76}
77
78
79template<class Thermo, class OtherThermo>
82(
83 const word& speciesName,
84 const volScalarField& Tf
85) const
86{
87 const label index = this->speciesNames_.find(speciesName);
88
89 if (index >= 0)
90 {
91 return
92 k_[index]
93 *this->otherThermo_.composition().Y(speciesName)
94 *this->otherThermo_.rhoThermo::rho()
95 /this->thermo_.rhoThermo::rho();
96 }
97 else
98 {
99 return
100 YSolvent_
101 *this->thermo_.composition().Y(speciesName);
102 }
103}
104
105
106template<class Thermo, class OtherThermo>
109(
110 const word& speciesName,
111 const volScalarField& Tf
112) const
113{
115 (
116 IOobject::groupName("YfPrime", this->pair_.name()),
117 this->pair_.phase1().mesh(),
119 );
120}
121
122
123// ************************************************************************* //
phaseModel & phase1
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Base class for interface composition models, templated on the two thermodynamic models either side of...
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Henry's law for gas solubility in liquid. The concentration of a dissolved species in the liquid is p...
Definition: Henry.H:64
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
Definition: Henry.C:82
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
Definition: Henry.C:109
virtual bool update()
Update the mesh for both mesh motion and topology change.
Description for mass transfer between a pair of phases. The direction of the mass transfer is from th...
Definition: phasePair.H:56
Lookup type of boundary radiation properties.
Definition: lookup.H:66
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
word timeName
Definition: getTimeIndex.H:3
const dimensionSet dimless
Dimensionless.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:54
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict