flowRateInletVelocityFvPatchVectorField.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) 2011-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 "one.H"
30 #include "patchPatchDist.H"
31 #include "wallPolyPatch.H"
32 
33 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
34 
35 void Foam::flowRateInletVelocityFvPatchVectorField::setWallDist()
36 {
37  if (profile_.valid())
38  {
39  const labelHashSet otherPatchIDs
40  (
41  patch().patch().boundaryMesh().findPatchIDs<wallPolyPatch>()
42  );
43 
44  const patchPatchDist pwd(patch().patch(), otherPatchIDs);
45 
46  y_ = pwd/gMax(pwd);
47  }
48 
49  area_ = gSum(patch().magSf());
50 }
51 
52 
54 Foam::flowRateInletVelocityFvPatchVectorField::profile()
55 {
56  if (profile_.valid())
57  {
58  return profile_->value(y_);
59  }
60  else
61  {
62  return tmp<scalarField>(new scalarField(size(), scalar(1)));
63  }
64 }
65 
66 
67 template<class ScaleType, class AlphaType, class RhoType>
68 void Foam::flowRateInletVelocityFvPatchVectorField::updateValues
69 (
70  const ScaleType& scale,
71  const AlphaType& alpha,
72  const RhoType& rho
73 )
74 {
75  const scalarField profile(this->profile());
76 
77  const scalar avgU =
78  -(scale*flowRate_->value(db().time().userTimeValue()))
79  /gSum(alpha*rho*profile*patch().magSf());
80 
81  operator==(avgU*profile*patch().nf());
82 }
83 
84 
85 template<class AlphaType>
86 void Foam::flowRateInletVelocityFvPatchVectorField::updateValues
87 (
88  const AlphaType& alpha
89 )
90 {
91  if (meanVelocity_)
92  {
93  updateValues(area_, alpha, one());
94  }
95  else if (volumetric_ || rhoName_ == "none")
96  {
97  updateValues(one(), alpha, one());
98  }
99  else
100  {
101  // Mass flow-rate
102  if (db().foundObject<volScalarField>(rhoName_))
103  {
104  const fvPatchField<scalar>& rhop =
105  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
106 
107  updateValues(one(), alpha, rhop);
108  }
109  else
110  {
111  // Use constant density
112  if (rhoInlet_ < 0)
113  {
115  << "Did not find registered density field " << rhoName_
116  << " and no constant density 'rhoInlet' specified"
117  << exit(FatalError);
118  }
119 
120  updateValues(one(), alpha, rhoInlet_);
121  }
122  }
123 }
124 
125 
126 bool Foam::flowRateInletVelocityFvPatchVectorField::canEvaluate()
127 {
128  return
130  || !patch().boundaryMesh().mesh().time().processorCase();
131 }
132 
133 
134 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
135 
138 (
139  const fvPatch& p,
141  const dictionary& dict
142 )
143 :
144  fixedValueFvPatchField<vector>(p, iF, dict, false),
145  rhoInlet_(dict.lookupOrDefault<scalar>("rhoInlet", -vGreat)),
146  alphaName_(dict.lookupOrDefault<word>("alpha", word::null))
147 {
148  if (dict.found("meanVelocity"))
149  {
150  meanVelocity_ = true;
151  volumetric_ = false;
152  flowRate_ = Function1<scalar>::New("meanVelocity", dict);
153  }
154  else if (dict.found("volumetricFlowRate"))
155  {
156  meanVelocity_ = false;
157  volumetric_ = true;
158  flowRate_ = Function1<scalar>::New("volumetricFlowRate", dict);
159  }
160  else if (dict.found("massFlowRate"))
161  {
162  meanVelocity_ = false;
163  volumetric_ = false;
164  flowRate_ = Function1<scalar>::New("massFlowRate", dict);
165  rhoName_ = word(dict.lookupOrDefault<word>("rho", "rho"));
166  }
167  else
168  {
170  (
171  dict
172  ) << "Please supply 'meanVelocity', 'volumetricFlowRate' or"
173  << " 'massFlowRate' and 'rho'" << exit(FatalIOError);
174  }
175 
176  if (dict.found("profile"))
177  {
178  profile_ = Function1<scalar>::New("profile", dict);
179  }
180 
181  if (canEvaluate())
182  {
183  setWallDist();
184  }
185 
186  if (!canEvaluate() || dict.found("value"))
187  {
189  (
190  vectorField("value", dict, p.size())
191  );
192  }
193  else
194  {
196  }
197 }
198 
199 
202 (
204  const fvPatch& p,
206  const fvPatchFieldMapper& mapper
207 )
208 :
209  fixedValueFvPatchField<vector>(ptf, p, iF, mapper),
210  flowRate_(ptf.flowRate_, false),
211  profile_(ptf.profile_, false),
212  meanVelocity_(ptf.meanVelocity_),
213  volumetric_(ptf.volumetric_),
214  rhoName_(ptf.rhoName_),
215  rhoInlet_(ptf.rhoInlet_),
216  alphaName_(ptf.alphaName_),
217  y_
218  (
219  profile_.valid() && canEvaluate()
220  ? mapper(ptf.y_)
221  : tmp<scalarField>(new scalarField())
222  ),
223  area_(NaN)
224 {}
225 
226 
229 (
232 )
233 :
234  fixedValueFvPatchField<vector>(ptf, iF),
235  flowRate_(ptf.flowRate_, false),
236  profile_(ptf.profile_, false),
237  meanVelocity_(ptf.meanVelocity_),
238  volumetric_(ptf.volumetric_),
239  rhoName_(ptf.rhoName_),
240  rhoInlet_(ptf.rhoInlet_),
241  alphaName_(ptf.alphaName_),
242  y_(ptf.y_),
243  area_(ptf.area_)
244 {}
245 
246 
247 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
248 
250 (
251  const fvPatchVectorField& ptf,
252  const fvPatchFieldMapper& mapper
253 )
254 {
255  fixedValueFvPatchVectorField::map(ptf, mapper);
256 
258  refCast<const flowRateInletVelocityFvPatchVectorField>(ptf);
259 
260  if (profile_.valid() && canEvaluate())
261  {
262  mapper(y_, tiptf.y_);
263  }
264 }
265 
266 
268 (
269  const fvPatchVectorField& ptf
270 )
271 {
272  fixedValueFvPatchVectorField::reset(ptf);
273 
275  refCast<const flowRateInletVelocityFvPatchVectorField>(ptf);
276 
277  if (profile_.valid() && canEvaluate())
278  {
279  y_.reset(tiptf.y_);
280  }
281 }
282 
283 
285 {
286  if (updated())
287  {
288  return;
289  }
290 
291  if (!canEvaluate())
292  {
294  << "Cannot evaluate flow rate on a non-parallel processor case"
295  << exit(FatalError);
296  }
297 
298  if (alphaName_ != word::null)
299  {
300  const fvPatchField<scalar>& alphap =
301  patch().lookupPatchField<volScalarField, scalar>(alphaName_);
302 
303  updateValues(alphap);
304  }
305  else
306  {
307  updateValues(one());
308  }
309 
310  fixedValueFvPatchVectorField::updateCoeffs();
311 }
312 
313 
315 {
317  writeEntry(os, flowRate_());
318  if (profile_.valid())
319  {
320  writeEntry(os, profile_());
321  }
322  if (!volumetric_)
323  {
324  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
325  writeEntryIfDifferent<scalar>(os, "rhoInlet", -vGreat, rhoInlet_);
326  }
327  writeEntryIfDifferent<word>(os, "alpha", word::null, alphaName_);
328  writeEntry(os, "value", *this);
329 }
330 
331 
332 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
333 
334 namespace Foam
335 {
337  (
340  );
341 }
342 
343 
344 // ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static autoPtr< Function1< Type > > New(const word &name, const dictionary &dict)
Selector.
Definition: Function1New.C:32
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
Velocity inlet boundary condition creating a velocity field with optionally specified profile normal ...
flowRateInletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchVectorField &, const fvPatchFieldMapper &)
Map the given fvPatchField onto this fvPatchField.
virtual void reset(const fvPatchVectorField &)
Reset the fvPatchField to the given fvPatchField.
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
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
A class representing the concept of 1 (scalar(1)) used to avoid unnecessary manipulations for objects...
Definition: one.H:51
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:318
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
bool valid(const PtrList< ModelType > &l)
static Type NaN()
Return a primitive with all components set to NaN.
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
Type gSum(const FieldField< Field, Type > &f)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:61
Field< vector > vectorField
Specialisation of Field<T> for vector.
IOerror FatalIOError
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:211
void evaluate(GeometricField< Type, PatchField, GeoMesh > &result, const Function1< Type > &func, const GeometricField< Type, PatchField, GeoMesh > &x)
error FatalError
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Type gMax(const FieldField< Field, Type > &f)
dictionary dict
volScalarField & p