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-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 
28 #include "volFields.H"
29 #include "surfaceFields.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  const dictionary& dict
80 )
81 :
82  mixedFvPatchScalarField(p, iF, dict, false),
83  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
84  UName_(dict.lookupOrDefault<word>("U", "U")),
85  phiYp_(p.size(), 0),
86  timeIndex_(-1),
87  c_(dict.lookupOrDefault<scalar>("c", scalar(0))),
88  property_
89  (
90  c_ == scalar(0)
91  ? massFraction
92  : propertyNames_.read(dict.lookup("property"))
93  )
94 {
96 
97  refValue() = Zero;
98  refGrad() = Zero;
99  valueFraction() = Zero;
100 }
101 
102 
105 (
107  const fvPatch& p,
109  const fvPatchFieldMapper& mapper
110 )
111 :
112  mixedFvPatchScalarField(ptf, p, iF, mapper),
113  phiName_(ptf.phiName_),
114  UName_(ptf.UName_),
115  phiYp_(p.size(), 0),
116  timeIndex_(-1),
117  c_(ptf.c_),
118  property_(ptf.property_)
119 {}
120 
121 
124 (
127 )
128 :
129  mixedFvPatchScalarField(ptf, iF),
130  phiName_(ptf.phiName_),
131  UName_(ptf.UName_),
132  phiYp_(ptf.size(), 0),
133  timeIndex_(-1),
134  c_(ptf.c_),
135  property_(ptf.property_)
136 {}
137 
138 
139 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
140 
142 (
143  const fvPatchScalarField& ptf,
144  const fvPatchFieldMapper& mapper
145 )
146 {
147  mixedFvPatchScalarField::map(ptf, mapper);
148 
150  refCast<const specieTransferMassFractionFvPatchScalarField>(ptf);
151 
152  mapper(phiYp_, tiptf.phiYp_);
153 }
154 
155 
157 (
158  const fvPatchScalarField& ptf
159 )
160 {
161  mixedFvPatchScalarField::reset(ptf);
162 
164  refCast<const specieTransferMassFractionFvPatchScalarField>(ptf);
165 
166  phiYp_.reset(tiptf.phiYp_);
167 }
168 
169 
170 const Foam::scalarField&
172 {
173  if (timeIndex_ != this->db().time().timeIndex())
174  {
175  timeIndex_ = this->db().time().timeIndex();
176 
177  phiYp_ = calcPhiYp();
178  }
179 
180  return phiYp_;
181 }
182 
183 
185 {
186  if (updated())
187  {
188  return;
189  }
190 
191  // Get the fluxes
192  const scalarField& phip =
193  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
194  const fvPatchVectorField& Up =
195  patch().lookupPatchField<volVectorField, vector>(UName_);
196  tmp<scalarField> uPhip =
197  refCast<const specieTransferVelocityFvPatchVectorField>(Up).phip();
198 
200  db().lookupType<fluidThermophysicalTransportModel>();
201 
202  // Get the diffusivity
203  const scalarField AAlphaEffp
204  (
205  patch().magSf()
206  *ttm.kappaEff(patch().index())
207  /ttm.thermo().Cp().boundaryField()[patch().index()]
208  );
209 
210  // Set the gradient and value so that the transport and diffusion combined
211  // result in the desired specie flux
212  valueFraction() = phip/(phip - patch().deltaCoeffs()*AAlphaEffp);
213  refValue() = *this;
214  refGrad() = phip*(*this - phiYp()/uPhip)/AAlphaEffp;
215 
216  mixedFvPatchScalarField::updateCoeffs();
217 }
218 
219 
221 (
222  Ostream& os
223 ) const
224 {
226  writeEntryIfDifferent<scalar>(os, "c", scalar(0), c_);
227  writeEntry(os, "property", propertyNames_[property_]);
228  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
229  writeEntryIfDifferent<word>(os, "U", "U", UName_);
230 }
231 
232 
233 // ************************************************************************* //
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Specialisation of basicMixture for a mixture consisting of a number for molecular species.
virtual const volScalarField & Cp() const =0
Heat capacity at constant pressure [J/kg/K].
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Base-class for multi-component fluid thermodynamic properties.
Abstract base class for fluid thermophysical transport models RAS, LES and laminar.
virtual const fluidThermo & thermo() const =0
Access function to fluid thermophysical properties.
Foam::fvPatchFieldMapper.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:87
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:231
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
Registry of regIOobjects.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
Abstract base class for specie-transferring mass fraction boundary conditions.
virtual const scalarField & phiYp() const
Return the flux of this species.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
static const basicSpecieMixture & composition(const objectRegistry &db)
Access the composition for the given database.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fvPatchFieldMapper &)
Map the given fvPatchField onto this fvPatchField.
static const NamedEnum< property, 4 > propertyNames_
Property type names.
specieTransferMassFractionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual tmp< volScalarField > kappaEff() const =0
Effective thermal turbulent conductivity.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
bool read(const char *, int32_t &)
Definition: int32IO.C:85
VolField< vector > volVectorField
Definition: volFieldsFwd.H:62
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
SurfaceField< scalar > surfaceScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
label timeIndex
Definition: getTimeIndex.H:4
dictionary dict
volScalarField & p
Foam::surfaceFields.