Saturated.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 | Copyright (C) 2015 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 "Saturated.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class Thermo, class OtherThermo>
33 wRatioByP() const
34 {
35  return
36  this->thermo_.composition().W(saturatedIndex_)
37  /this->thermo_.composition().W()
38  /this->thermo_.p();
39 }
40 
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
44 template<class Thermo, class OtherThermo>
46 (
47  const dictionary& dict,
48  const phasePair& pair
49 )
50 :
51  InterfaceCompositionModel<Thermo, OtherThermo>(dict, pair),
52  saturatedName_(this->speciesNames_[0]),
53  saturatedIndex_
54  (
55  this->thermo_.composition().species()[saturatedName_]
56  ),
57  saturationModel_
58  (
59  saturationModel::New
60  (
61  dict.subDict("saturationPressure")
62  )
63  )
64 {
65  if (this->speciesNames_.size() != 1)
66  {
68  (
69  "template<class Thermo, class OtherThermo>"
70  "Foam::interfaceCompositionModels::Saturated<Thermo, OtherThermo>::"
71  "Saturated"
72  "( "
73  "const dictionary& dict, "
74  "const phasePair& pair "
75  ")"
76  ) << "Saturated model is suitable for one species only."
77  << exit(FatalError);
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
83 
84 template<class Thermo, class OtherThermo>
86 {}
87 
88 
89 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
90 
91 template<class Thermo, class OtherThermo>
92 void
94 (
95  const volScalarField& Tf
96 )
97 {
98  // do nothing
99 }
100 
101 
102 template<class Thermo, class OtherThermo>
105 (
106  const word& speciesName,
107  const volScalarField& Tf
108 ) const
109 {
110  if (saturatedName_ == speciesName)
111  {
112  return wRatioByP()*saturationModel_->pSat(Tf);
113  }
114  else
115  {
116  const label speciesIndex
117  (
118  this->thermo_.composition().species()[speciesName]
119  );
120 
121  return
122  this->thermo_.Y()[speciesIndex]
123  *(scalar(1) - wRatioByP()*saturationModel_->pSat(Tf))
124  /max(scalar(1) - this->thermo_.Y()[saturatedIndex_], SMALL);
125  }
126 }
127 
128 
129 template<class Thermo, class OtherThermo>
132 (
133  const word& speciesName,
134  const volScalarField& Tf
135 ) const
136 {
137  if (saturatedName_ == speciesName)
138  {
139  return wRatioByP()*saturationModel_->pSatPrime(Tf);
140  }
141  else
142  {
143  const label speciesIndex
144  (
145  this->thermo_.composition().species()[speciesName]
146  );
147 
148  return
149  - this->thermo_.Y()[speciesIndex]
150  *wRatioByP()*saturationModel_->pSatPrime(Tf)
151  /max(scalar(1) - this->thermo_.Y()[saturatedIndex_], SMALL);
152  }
153 }
154 
155 
156 // ************************************************************************* //
virtual tmp< volScalarField > Yf(const word &speciesName, const volScalarField &Tf) const
The interface species fraction.
Saturated(const dictionary &dict, const phasePair &pair)
Construct from components.
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
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual tmp< volScalarField > YfPrime(const word &speciesName, const volScalarField &Tf) const
The interface species fraction derivative w.r.t. temperature.
tmp< volScalarField > wRatioByP() const
Constant of propotionality between partial pressure and mass.
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
basicMultiComponentMixture & composition
Definition: createFields.H:35
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
const hashedWordList speciesNames_
Names of the transferring species.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
virtual void update(const volScalarField &Tf)
Update the composition.
error FatalError
A class for managing temporary objects.
Definition: PtrList.H:118