totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9  Copyright (C) 2018-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
31 #include "fvPatchFieldMapper.H"
32 #include "volFields.H"
33 #include "surfaceFields.H"
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
40 (
41  const fvPatch& p,
43 )
44 :
46  phiName_("phi"),
47  rhoName_("none"),
48  massFluxFraction_(1.0)
49 {
50  refValue() = 0.0;
51  refGrad() = 0.0;
52  valueFraction() = 0.0;
53 }
54 
55 
58 (
59  const fvPatch& p,
61  const dictionary& dict
62 )
63 :
65  phiName_(dict.getOrDefault<word>("phi", "phi")),
66  rhoName_(dict.getOrDefault<word>("rho", "none")),
67  massFluxFraction_(dict.getOrDefault<scalar>("massFluxFraction", 1))
68 {
69 
70  refValue() = 1.0;
71  refGrad() = 0.0;
72  valueFraction() = 0.0;
73 
74  if (dict.found("value"))
75  {
77  (
78  Field<scalar>("value", dict, p.size())
79  );
80  }
81  else
82  {
84  }
85 }
86 
87 
90 (
92  const fvPatch& p,
94  const fvPatchFieldMapper& mapper
95 )
96 :
97  mixedFvPatchField<scalar>(ptf, p, iF, mapper),
98  phiName_(ptf.phiName_),
99  rhoName_(ptf.rhoName_),
100  massFluxFraction_(ptf.massFluxFraction_)
101 {}
102 
103 
106 (
108 )
109 :
111  phiName_(tppsf.phiName_),
112  rhoName_(tppsf.rhoName_),
113  massFluxFraction_(tppsf.massFluxFraction_)
114 {}
115 
116 
119 (
122 )
123 :
124  mixedFvPatchField<scalar>(tppsf, iF),
125  phiName_(tppsf.phiName_),
126  rhoName_(tppsf.rhoName_),
127  massFluxFraction_(tppsf.massFluxFraction_)
128 {}
129 
130 
131 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
132 
134 (
135  const fvPatchFieldMapper& m
136 )
137 {
138  scalarField::autoMap(m);
139 }
140 
141 
143 (
144  const fvPatchScalarField& ptf,
145  const labelList& addr
146 )
147 {
149 }
150 
151 
153 {
154  if (this->updated())
155  {
156  return;
157  }
158 
159  const label patchi = patch().index();
160 
161  const compressible::turbulenceModel& turbModel =
163  (
165  (
167  internalField().group()
168  )
169  );
170 
171  const fvsPatchField<scalar>& phip =
172  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
173 
174  const scalarField alphap(turbModel.alphaEff(patchi));
175 
176  refValue() = massFluxFraction_;
177  refGrad() = 0.0;
178 
179  valueFraction() =
180  1.0
181  /(
182  1.0
183  + alphap*patch().deltaCoeffs()*patch().magSf()/max(mag(phip), SMALL)
184  );
185 
187 
188  if (debug)
189  {
190  scalar phi = gSum(-phip*(*this));
191 
192  Info<< patch().boundaryMesh().mesh().name() << ':'
193  << patch().name() << ':'
194  << this->internalField().name() << " :"
195  << " mass flux[Kg/s]:" << phi
196  << endl;
197  }
198 }
199 
200 
202 (
203  Ostream& os
204 ) const
205 {
207  os.writeEntry("phi", phiName_);
208  os.writeEntry("rho", rhoName_);
209  os.writeEntry("massFluxFraction", massFluxFraction_);
210  this->writeEntry("value", os);
211 }
212 
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 namespace Foam
217 {
219  (
222  );
223 
224 }
225 
226 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::fvPatchField
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: volSurfaceMapping.H:51
volFields.H
totalFlowRateAdvectiveDiffusiveFvPatchScalarField.H
Foam::fvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::fvPatchField< scalar >::internalField
const DimensionedField< scalar, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:363
Foam::mixedFvPatchField< scalar >::valueFraction
virtual scalarField & valueFraction()
Definition: mixedFvPatchField.H:250
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::constant::atomic::group
constexpr const char *const group
Group name for atomic constants.
Definition: atomicConstants.H:52
Foam::mixedFvPatchField< scalar >::refValue
virtual Field< scalar > & refValue()
Definition: mixedFvPatchField.H:230
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField::rmap
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C:143
Foam::fvsPatchField
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:68
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField::totalFlowRateAdvectiveDiffusiveFvPatchScalarField
totalFlowRateAdvectiveDiffusiveFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C:40
Foam::mixedFvPatchField< scalar >::refGrad
virtual Field< scalar > & refGrad()
Definition: mixedFvPatchField.H:240
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
surfaceFields.H
Foam::surfaceFields.
fvPatchFieldMapper.H
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::fvPatch::name
virtual const word & name() const
Return name.
Definition: fvPatch.H:167
Foam::fvBoundaryMesh::mesh
const fvMesh & mesh() const noexcept
Return the mesh reference.
Definition: fvBoundaryMesh.H:103
Foam::fvPatchField< scalar >::updated
bool updated() const
Return true if the boundary condition has already been updated.
Definition: fvPatchField.H:387
Foam::fvPatch::boundaryMesh
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:203
Foam::Field< scalar >
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C:202
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::objectRegistry::lookupObject
const Type & lookupObject(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:434
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::fvPatch::lookupPatchField
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
Definition: fvPatchFvMeshTemplates.C:34
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::fvPatch::index
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:197
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField
This BC is used for species inlets. The diffusion and advection fluxes are considered to calculate th...
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.H:53
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::mixedFvPatchField< scalar >
Foam::fvPatchField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:321
Foam::ThermalDiffusivity::alphaEff
virtual tmp< volScalarField > alphaEff() const
Return the effective turbulent thermal diffusivity for enthalpy.
Definition: ThermalDiffusivity.H:160
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition: surfaceFieldsFwd.H:54
Foam::fvPatchField< scalar >::db
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:206
Foam::ThermalDiffusivity
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
Definition: phaseCompressibleTurbulenceModelFwd.H:47
Foam::List< label >
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::fvPatch::deltaCoeffs
const scalarField & deltaCoeffs() const
Definition: fvPatch.C:196
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField::autoMap
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C:134
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::fvPatchField< scalar >::patch
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:357
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::fvPatch::magSf
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:156
Foam::totalFlowRateAdvectiveDiffusiveFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: totalFlowRateAdvectiveDiffusiveFvPatchScalarField.C:152
turbulentFluidThermoModel.H
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::fvMesh::name
const word & name() const
Return reference to name.
Definition: fvMesh.H:300