specieTransferMassFractionFvPatchScalarField.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) 2019-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 
28 #include "volFields.H"
29 #include "surfaceFields.H"
31 #include "psiReactionThermo.H"
32 #include "rhoReactionThermo.H"
33 
34 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  template<>
39  const char* NamedEnum
40  <
42  4
43  >::names[] =
44  {
45  "massFraction",
46  "moleFraction",
47  "molarConcentration",
48  "partialPressure"
49  };
50 }
51 
52 const Foam::NamedEnum
53 <
55  4
57 
58 
59 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
60 
63 (
64  const objectRegistry& db
65 )
66 {
68 
69  if (db.foundObject<psiReactionThermo>(name))
70  {
72  }
73  else if (db.foundObject<rhoReactionThermo>(name))
74  {
76  }
77  else
78  {
80  << "Could not find a multi-component thermodynamic model."
81  << exit(FatalError);
82 
83  return NullObjectRef<basicSpecieMixture>();
84  }
85 }
86 
87 
88 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
89 
92 (
93  const fvPatch& p,
95 )
96 :
97  mixedFvPatchScalarField(p, iF),
98  phiName_("phi"),
99  UName_("U"),
100  phiYp_(p.size(), 0),
101  timeIndex_(-1),
102  c_(0),
103  property_(massFraction)
104 {
105  refValue() = Zero;
106  refGrad() = Zero;
107  valueFraction() = Zero;
108 }
109 
110 
113 (
114  const fvPatch& p,
116  const dictionary& dict
117 )
118 :
119  mixedFvPatchScalarField(p, iF),
120  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
121  UName_(dict.lookupOrDefault<word>("U", "U")),
122  phiYp_(p.size(), 0),
123  timeIndex_(-1),
124  c_(dict.lookupOrDefault<scalar>("c", scalar(0))),
125  property_
126  (
127  c_ == scalar(0)
128  ? massFraction
129  : propertyNames_.read(dict.lookup("property"))
130  )
131 {
132  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
133 
134  refValue() = Zero;
135  refGrad() = Zero;
136  valueFraction() = Zero;
137 }
138 
139 
142 (
144  const fvPatch& p,
146  const fvPatchFieldMapper& mapper
147 )
148 :
149  mixedFvPatchScalarField(ptf, p, iF, mapper),
150  phiName_(ptf.phiName_),
151  UName_(ptf.UName_),
152  phiYp_(p.size(), 0),
153  timeIndex_(-1),
154  c_(ptf.c_),
155  property_(ptf.property_)
156 {}
157 
158 
161 (
163 )
164 :
165  mixedFvPatchScalarField(ptf),
166  phiName_(ptf.phiName_),
167  UName_(ptf.UName_),
168  phiYp_(ptf.size(), 0),
169  timeIndex_(-1),
170  c_(ptf.c_),
171  property_(ptf.property_)
172 {}
173 
174 
177 (
180 )
181 :
182  mixedFvPatchScalarField(ptf, iF),
183  phiName_(ptf.phiName_),
184  UName_(ptf.UName_),
185  phiYp_(ptf.size(), 0),
186  timeIndex_(-1),
187  c_(ptf.c_),
188  property_(ptf.property_)
189 {}
190 
191 
192 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
193 
195 (
196  const fvPatchFieldMapper& m
197 )
198 {
199  mixedFvPatchScalarField::autoMap(m);
200 
201  m(phiYp_, phiYp_);
202 }
203 
204 
206 (
207  const fvPatchScalarField& ptf,
208  const labelList& addr
209 )
210 {
211  mixedFvPatchScalarField::rmap(ptf, addr);
212 
214  refCast<const specieTransferMassFractionFvPatchScalarField>(ptf);
215 
216  phiYp_.rmap(tiptf.phiYp_, addr);
217 }
218 
219 
220 const Foam::scalarField&
222 {
223  if (timeIndex_ != this->db().time().timeIndex())
224  {
225  timeIndex_ = this->db().time().timeIndex();
226 
227  phiYp_ = calcPhiYp();
228  }
229 
230  return phiYp_;
231 }
232 
233 
235 {
236  if (updated())
237  {
238  return;
239  }
240 
241  // Get the fluxes
242  const scalarField& phip =
243  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
244  const fvPatchVectorField& Up =
245  patch().lookupPatchField<volVectorField, vector>(UName_);
246  tmp<scalarField> uPhip =
247  refCast<const specieTransferVelocityFvPatchVectorField>(Up).phip();
248 
249  // Get the diffusivity
250  // !!! <-- This is a potential lagging issue as alphaEff(patchi) calculates
251  // alpha on demand, so the value may be different from alphaEff() which
252  // uses the alpha field cached by the thermodynamics.
253  const scalarField AAlphaEffp
254  (
255  patch().magSf()
256  *db().lookupObject<thermophysicalTransportModel>
257  (
258  thermophysicalTransportModel::typeName
259  )
260  .alphaEff(patch().index())
261  );
262 
263  // Set the gradient and value so that the transport and diffusion combined
264  // result in the desired specie flux
265  valueFraction() = phip/(phip - patch().deltaCoeffs()*AAlphaEffp);
266  refValue() = *this;
267  refGrad() = phip*(*this - phiYp()/uPhip)/AAlphaEffp;
268 
269  mixedFvPatchScalarField::updateCoeffs();
270 }
271 
272 
274 (
275  Ostream& os
276 ) const
277 {
279  writeEntryIfDifferent<scalar>(os, "c", scalar(0), c_);
280  writeEntry(os, "property", propertyNames_[property_]);
281  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
282  writeEntryIfDifferent<word>(os, "U", "U", UName_);
283 }
284 
285 
286 // ************************************************************************* //
Foam::surfaceFields.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
basicSpecieMixture & composition
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
bool foundObject(const word &name) const
Is the named Type found?
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Specialization of basicMixture for a mixture consisting of a number for molecular species...
static const basicSpecieMixture & composition(const objectRegistry &db)
Access the composition for the given database.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
specieTransferMassFractionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::rhoReactionThermo.
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:97
virtual label size() const
Return size.
Definition: fvPatch.H:155
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
virtual const scalarField & phiYp() const
Return the flux of this species.
Foam::psiReactionThermo.
const Time & time() const
Return time.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
Abstract base class for specie-transferring mass fraction boundary conditions.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:295
label timeIndex
Definition: getTimeIndex.H:4
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Registry of regIOobjects.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
static const NamedEnum< property, 4 > propertyNames_
Property type names.
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:812