semiPermeableBaffleMassFractionFvPatchScalarField.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) 2017-2019 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 "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
31 #include "turbulenceModel.H"
32 #include "psiReactionThermo.H"
33 #include "rhoReactionThermo.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  template<>
40  const char* NamedEnum
41  <
43  4
44  >::names[] =
45  {
46  "none",
47  "massFraction",
48  "moleFraction",
49  "partialPressure",
50  };
51 }
52 
53 const Foam::NamedEnum
54 <
56  4
58 
59 
60 // * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
61 
64 (
65  const objectRegistry& db
66 )
67 {
69 
70  if (db.foundObject<psiReactionThermo>(name))
71  {
73  }
74  else if (db.foundObject<rhoReactionThermo>(name))
75  {
77  }
78  else
79  {
81  << "Could not find a multi-component thermodynamic model."
82  << exit(FatalError);
83 
84  return NullObjectRef<basicSpecieMixture>();
85  }
86 }
87 
88 
89 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
90 
93 (
94  const fvPatch& p,
96 )
97 :
99  mixedFvPatchScalarField(p, iF),
100  c_(0),
101  input_(none),
102  phiName_("phi"),
103  pName_("p")
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  mappedPatchBase(p.patch(), NEARESTPATCHFACE, dict),
120  mixedFvPatchScalarField(p, iF),
121  c_(dict.lookupOrDefault<scalar>("c", scalar(0))),
122  input_
123  (
124  c_ == scalar(0)
125  ? none
126  : inputNames_.read(dict.lookup("input"))
127  ),
128  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
129  pName_(dict.lookupOrDefault<word>("p", "p"))
130 {
131  fvPatchScalarField::operator=(scalarField("value", dict, p.size()));
132 
133  refValue() = Zero;
134  refGrad() = Zero;
135  valueFraction() = Zero;
136 }
137 
138 
141 (
143  const fvPatch& p,
145  const fvPatchFieldMapper& mapper
146 )
147 :
148  mappedPatchBase(p.patch(), ptf),
149  mixedFvPatchScalarField(ptf, p, iF, mapper),
150  c_(ptf.c_),
151  input_(ptf.input_),
152  phiName_(ptf.phiName_),
153  pName_(ptf.pName_)
154 {}
155 
156 
159 (
161 )
162 :
163  mappedPatchBase(ptf.patch().patch(), ptf),
164  mixedFvPatchScalarField(ptf),
165  c_(ptf.c_),
166  input_(ptf.input_),
167  phiName_(ptf.phiName_),
168  pName_(ptf.pName_)
169 {}
170 
171 
174 (
177 )
178 :
179  mappedPatchBase(ptf.patch().patch(), ptf),
180  mixedFvPatchScalarField(ptf, iF),
181  c_(ptf.c_),
182  input_(ptf.input_),
183  phiName_(ptf.phiName_),
184  pName_(ptf.pName_)
185 {}
186 
187 
188 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
189 
192 {
193  if (c_ == scalar(0))
194  {
195  return tmp<scalarField>(new scalarField(patch().size(), Zero));
196  }
197 
198  const word& YName = internalField().name();
199 
200  // Initialise the input variables to the mass fractions
201  scalarField psic(patchInternalField());
202 
203  const label nbrPatchi = samplePolyPatch().index();
204  const fvPatch& nbrPatch = patch().boundaryMesh()[nbrPatchi];
205  const fvPatchScalarField& nbrYp =
206  nbrPatch.lookupPatchField<volScalarField, scalar>(YName);
207  scalarField nbrPsic(nbrYp.patchInternalField());
208  mappedPatchBase::map().distribute(nbrPsic);
209 
210  switch (input_)
211  {
212  case none:
214  << "A none input cannot be used with a non-zero transfer "
215  << "coefficient" << exit(FatalError);
216 
217  case massFraction:
218  // Do nothing
219  break;
220 
221  case partialPressure:
222  // Multiply by pressure
223  {
224  psic *=
225  patch().lookupPatchField<volScalarField, scalar>(pName_);
226 
227  fvPatchScalarField nbrP
228  (
229  nbrPatch.lookupPatchField<volScalarField, scalar>(pName_)
230  );
232  nbrPsic *= nbrP;
233  }
234 
235  // Falls through ...
236 
237  case moleFraction:
238  // Convert to mole fraction
239  {
240  const basicSpecieMixture& mixture = composition(db());
241  const scalar Wi(mixture.Wi(mixture.species()[YName]));
242  const basicThermo& thermo =
244 
245  psic *= thermo.W(patch().index())/Wi;
246 
247  scalarField nbrW(thermo.W(nbrPatch.index()));
249  nbrPsic *= nbrW/Wi;
250  }
251  break;
252  }
253 
254  return c_*patch().magSf()*(psic - nbrPsic);
255 }
256 
257 
259 {
260  if (updated())
261  {
262  return;
263  }
264 
265  const scalarField& phip =
266  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
267 
268  const turbulenceModel& turbModel =
270  (
272  );
273  const scalarField muEffp(turbModel.muEff(patch().index()));
274  const scalarField AMuEffp(patch().magSf()*muEffp);
275 
276  valueFraction() = phip/(phip - patch().deltaCoeffs()*AMuEffp);
277  refGrad() = - phiY()/AMuEffp;
278 
279  mixedFvPatchScalarField::updateCoeffs();
280 }
281 
282 
284 (
285  Ostream& os
286 ) const
287 {
290  writeEntryIfDifferent<scalar>(os, "c", scalar(0), c_);
291  writeEntry(os, "input", inputNames_[input_]);
292  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
293  writeEntryIfDifferent<word>(os, "p", "p", pName_);
294  writeEntry(os, "value", *this);
295 }
296 
297 
298 // * * * * * * * * * * * * * * Build Macro Function * * * * * * * * * * * * //
299 
300 namespace Foam
301 {
303  (
306  );
307 }
308 
309 // ************************************************************************* //
Foam::surfaceFields.
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:179
dictionary dict
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
Abstract base-class for fluid and solid thermodynamic properties.
Definition: basicThermo.H:52
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:173
semiPermeableBaffleMassFractionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
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
void distribute(List< T > &fld, const bool dummyTransform=true, const int tag=UPstream::msgType()) const
Distribute data using default commsType.
virtual void write(Ostream &) const
Write as a dictionary.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
const speciesTable & species() const
Return the table of species.
This is a mass-fraction boundary condition for a semi-permeable baffle.
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 scalar Wi(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
rhoReactionThermo & thermo
Definition: createFields.H:28
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Specialization of basicMultiComponentMixture for a mixture consisting of a number for molecular speci...
const mapDistribute & map() const
Return reference to the parallel distribution map.
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
Abstract base class for turbulence models (RAS, LES and laminar).
Macros for easy insertion into run-time selection tables.
const word dictName() const
Return the local dictionary name (final part of scoped name)
Definition: dictionary.H:123
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Foam::rhoReactionThermo.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:137
Foam::fvPatchFieldMapper.
Determines a mapping between patch face centres and mesh cell or face centres and processors they&#39;re ...
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:53
Foam::psiReactionThermo.
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:194
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
tmp< scalarField > phiY() const
Return the flux of this species through the baffle.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
static const NamedEnum< input, 4 > inputNames_
Input variable type names.
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
Lookup and return the patchField of the named field from the.
static const basicSpecieMixture & composition(const objectRegistry &db)
Access the composition for the given database.
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:295
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual tmp< volScalarField > W() const =0
Molecular weight [kg/kmol].
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Registry of regIOobjects.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583