specieTransferTemperatureFvPatchScalarField.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 
29 #include "volFields.H"
30 #include "surfaceFields.H"
32 #include "basicSpecieMixture.H"
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
39 (
40  const fvPatch& p,
42 )
43 :
45  phiName_("phi"),
46  UName_("U")
47 {}
48 
49 
52 (
53  const fvPatch& p,
55  const dictionary& dict,
56  const bool readValue
57 )
58 :
60  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
61  UName_(dict.lookupOrDefault<word>("U", "U"))
62 {
63  if (readValue)
64  {
65  fvPatchScalarField::operator==(scalarField("value", dict, p.size()));
66  }
67 }
68 
69 
72 (
74  const fvPatch& p,
76  const fvPatchFieldMapper& mapper
77 )
78 :
80  phiName_(ptf.phiName_),
81  UName_(ptf.UName_)
82 {}
83 
84 
87 (
90 )
91 :
93  phiName_(ptf.phiName_),
94  UName_(ptf.UName_)
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
102 {
105  const PtrList<volScalarField>& Y = mixture.Y();
106 
107  // Get thermodynamic properties
108  const fluidThermo& thermo =
109  db().lookupObject<fluidThermo>(basicThermo::dictName);
110  const fvPatchScalarField& Tp = *this;
111  const fvPatchScalarField& pp = thermo.p().boundaryField()[patch().index()];
112 
113  // Sum up the phiHep from all the species
114  tmp<scalarField> tPhiHep(new scalarField(this->size(), 0));
115  scalarField& phiHep = tPhiHep.ref();
116  forAll(Y, i)
117  {
118  const fvPatchScalarField& Yp = Y[i].boundaryField()[patch().index()];
119 
120  if (!isA<YBCType>(Yp))
121  {
123  << "The mass-fraction condition on patch " << patch().name()
124  << " is not of type " << YBCType::typeName << "."
125  << exit(FatalError);
126  }
127 
128  phiHep += refCast<const YBCType>(Yp).phiYp()*mixture.HE(i, pp, Tp);
129  }
130 
131  return tPhiHep;
132 }
133 
134 
136 {
137  if (updated())
138  {
139  return;
140  }
141 
142  // Get the fluxes
143  const scalarField& phip =
144  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
145  const fvPatchVectorField& Up =
146  patch().lookupPatchField<volVectorField, vector>(UName_);
147  tmp<scalarField> uPhip =
148  refCast<const specieTransferVelocityFvPatchVectorField>(Up).phip();
149 
150  // Get the diffusivity
151  // !!! <-- This is a potential lagging issue as alphaEff(patchi) calculates
152  // alpha on demand, so the value may be different from alphaEff() which
153  // uses the alpha field cached by the thermodynamics
154  const scalarField AAlphaEffp
155  (
156  patch().magSf()
157  *db().lookupObject<thermophysicalTransportModel>
158  (
159  thermophysicalTransportModel::typeName
160  ).alphaEff(patch().index())
161  );
162 
163  // Get the current energy to linearise around
164  const fluidThermo& thermo =
165  db().lookupObject<fluidThermo>(basicThermo::dictName);
166  const scalarField& hep = thermo.he().boundaryField()[patch().index()];
167 
168  heValueFraction() = phip/(phip - patch().deltaCoeffs()*AAlphaEffp);
169  heRefValue() = hep;
170  heRefGrad() = phip*(hep - phiHep()/uPhip)/AAlphaEffp;
171 
172  mixedEnergyCalculatedTemperatureFvPatchScalarField::updateCoeffs();
173 }
174 
175 
177 (
178  Ostream& os
179 ) const
180 {
182  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
183  writeEntryIfDifferent<word>(os, "U", "U", UName_);
184  writeEntry(os, "value", *this);
185 }
186 
187 
188 // * * * * * * * * * * * * * * Build Macro Function * * * * * * * * * * * * //
189 
190 namespace Foam
191 {
193  (
196  );
197 }
198 
199 // ************************************************************************* //
Foam::surfaceFields.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
fluidReactionThermo & thermo
Definition: createFields.H:28
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:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
basicSpecieMixture & composition
const Boundary & boundaryField() const
Return const-reference to the boundary field.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
specieTransferTemperatureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
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:260
Specialisation of basicMixture for a mixture consisting of a number for molecular species...
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:58
virtual volScalarField & p()=0
Pressure [Pa].
Macros for easy insertion into run-time selection tables.
virtual volScalarField & he()=0
Enthalpy/Internal energy [J/kg].
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
Foam::fvPatchFieldMapper.
virtual scalar HE(const label speciei, const scalar p, const scalar T) const =0
Enthalpy/Internal energy [J/kg].
virtual label size() const
Return size.
Definition: fvPatch.H:156
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
volScalarField scalarField(fieldObject, mesh)
static const word dictName
Name of the thermophysical properties dictionary.
Definition: basicThermo.H:129
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
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,.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
Base class for temperature boundary conditions in which the parameters of the mixed energy condition ...
PtrList< volScalarField > & Y
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
This is a temperature boundary condition for a specie-transferring wall.
const tmp< scalarField > phiHep() const
Return the flux of energy.
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Namespace for OpenFOAM.