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-2022 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 "fluidReactionThermo.H"
32 
33 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  template<>
38  const char* NamedEnum
39  <
41  4
42  >::names[] =
43  {
44  "massFraction",
45  "moleFraction",
46  "molarConcentration",
47  "partialPressure"
48  };
49 }
50 
51 const Foam::NamedEnum
52 <
54  4
56 
57 
58 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
59 
62 (
63  const objectRegistry& db
64 )
65 {
66  const word& name = physicalProperties::typeName;
67 
69 }
70 
71 
72 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
73 
76 (
77  const fvPatch& p,
79 )
80 :
81  mixedFvPatchScalarField(p, iF),
82  phiName_("phi"),
83  UName_("U"),
84  phiYp_(p.size(), 0),
85  timeIndex_(-1),
86  c_(0),
87  property_(massFraction)
88 {
89  refValue() = Zero;
90  refGrad() = Zero;
91  valueFraction() = Zero;
92 }
93 
94 
97 (
98  const fvPatch& p,
100  const dictionary& dict
101 )
102 :
103  mixedFvPatchScalarField(p, iF),
104  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
105  UName_(dict.lookupOrDefault<word>("U", "U")),
106  phiYp_(p.size(), 0),
107  timeIndex_(-1),
108  c_(dict.lookupOrDefault<scalar>("c", scalar(0))),
109  property_
110  (
111  c_ == scalar(0)
112  ? massFraction
113  : propertyNames_.read(dict.lookup("property"))
114  )
115 {
116  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
117 
118  refValue() = Zero;
119  refGrad() = Zero;
120  valueFraction() = Zero;
121 }
122 
123 
126 (
128  const fvPatch& p,
130  const fvPatchFieldMapper& mapper
131 )
132 :
133  mixedFvPatchScalarField(ptf, p, iF, mapper),
134  phiName_(ptf.phiName_),
135  UName_(ptf.UName_),
136  phiYp_(p.size(), 0),
137  timeIndex_(-1),
138  c_(ptf.c_),
139  property_(ptf.property_)
140 {}
141 
142 
145 (
148 )
149 :
150  mixedFvPatchScalarField(ptf, iF),
151  phiName_(ptf.phiName_),
152  UName_(ptf.UName_),
153  phiYp_(ptf.size(), 0),
154  timeIndex_(-1),
155  c_(ptf.c_),
156  property_(ptf.property_)
157 {}
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
163 (
164  const fvPatchFieldMapper& m
165 )
166 {
167  mixedFvPatchScalarField::autoMap(m);
168 
169  m(phiYp_, phiYp_);
170 }
171 
172 
174 (
175  const fvPatchScalarField& ptf,
176  const labelList& addr
177 )
178 {
179  mixedFvPatchScalarField::rmap(ptf, addr);
180 
182  refCast<const specieTransferMassFractionFvPatchScalarField>(ptf);
183 
184  phiYp_.rmap(tiptf.phiYp_, addr);
185 }
186 
187 
189 (
190  const fvPatchScalarField& ptf
191 )
192 {
193  mixedFvPatchScalarField::reset(ptf);
194 
196  refCast<const specieTransferMassFractionFvPatchScalarField>(ptf);
197 
198  phiYp_.reset(tiptf.phiYp_);
199 }
200 
201 
202 const Foam::scalarField&
204 {
205  if (timeIndex_ != this->db().time().timeIndex())
206  {
207  timeIndex_ = this->db().time().timeIndex();
208 
209  phiYp_ = calcPhiYp();
210  }
211 
212  return phiYp_;
213 }
214 
215 
217 {
218  if (updated())
219  {
220  return;
221  }
222 
223  // Get the fluxes
224  const scalarField& phip =
225  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
226  const fvPatchVectorField& Up =
227  patch().lookupPatchField<volVectorField, vector>(UName_);
228  tmp<scalarField> uPhip =
229  refCast<const specieTransferVelocityFvPatchVectorField>(Up).phip();
230 
231  // Get the diffusivity
232  // !!! <-- This is a potential lagging issue as alphaEff(patchi) calculates
233  // alpha on demand, so the value may be different from alphaEff() which
234  // uses the alpha field cached by the thermodynamics.
235  const scalarField AAlphaEffp
236  (
237  patch().magSf()
238  *db().lookupObject<thermophysicalTransportModel>
239  (
240  thermophysicalTransportModel::typeName
241  )
242  .alphaEff(patch().index())
243  );
244 
245  // Set the gradient and value so that the transport and diffusion combined
246  // result in the desired specie flux
247  valueFraction() = phip/(phip - patch().deltaCoeffs()*AAlphaEffp);
248  refValue() = *this;
249  refGrad() = phip*(*this - phiYp()/uPhip)/AAlphaEffp;
250 
251  mixedFvPatchScalarField::updateCoeffs();
252 }
253 
254 
256 (
257  Ostream& os
258 ) const
259 {
261  writeEntryIfDifferent<scalar>(os, "c", scalar(0), c_);
262  writeEntry(os, "property", propertyNames_[property_]);
263  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
264  writeEntryIfDifferent<word>(os, "U", "U", UName_);
265 }
266 
267 
268 // ************************************************************************* //
Foam::surfaceFields.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
basicSpecieMixture & composition
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:63
Base-class for multi-component fluid thermodynamic properties.
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:243
Specialisation 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:59
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
specieTransferMassFractionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:97
virtual label size() const
Return size.
Definition: fvPatch.H:157
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.
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:28
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:258
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:864