flowRateInletVelocityFvPatchVectorField.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-2017 OpenFOAM Foundation
9  Copyright (C) 2019-2021 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 "volFields.H"
32 #include "one.H"
33 #include "Switch.H"
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
39 (
40  const fvPatch& p,
42 )
43 :
45  flowRate_(),
46  rhoName_("rho"),
47  rhoInlet_(0.0),
48  volumetric_(false),
49  extrapolateProfile_(false)
50 {}
51 
52 
55 (
56  const fvPatch& p,
58  const dictionary& dict
59 )
60 :
62  rhoName_("rho"),
63  rhoInlet_(dict.getOrDefault<scalar>("rhoInlet", -VGREAT)),
64  volumetric_(false),
65  extrapolateProfile_
66  (
67  dict.getOrDefault<Switch>("extrapolateProfile", false)
68  )
69 {
70  if (dict.found("volumetricFlowRate"))
71  {
72  volumetric_ = true;
73  flowRate_ =
74  Function1<scalar>::New("volumetricFlowRate", dict, &db());
75  }
76  else if (dict.found("massFlowRate"))
77  {
78  volumetric_ = false;
79  flowRate_ = Function1<scalar>::New("massFlowRate", dict, &db());
80  rhoName_ = dict.getOrDefault<word>("rho", "rho");
81  }
82  else
83  {
85  << "Please supply either 'volumetricFlowRate' or"
86  << " 'massFlowRate' and 'rho'" << nl
87  << exit(FatalIOError);
88  }
89 
90  // Value field require if mass based
91  if (dict.found("value"))
92  {
94  (
95  vectorField("value", dict, p.size())
96  );
97  }
98  else
99  {
100  evaluate(Pstream::commsTypes::blocking);
101  }
102 }
103 
104 
107 (
109  const fvPatch& p,
111  const fvPatchFieldMapper& mapper
112 )
113 :
114  fixedValueFvPatchField<vector>(ptf, p, iF, mapper),
115  flowRate_(ptf.flowRate_.clone()),
116  rhoName_(ptf.rhoName_),
117  rhoInlet_(ptf.rhoInlet_),
118  volumetric_(ptf.volumetric_),
119  extrapolateProfile_(ptf.extrapolateProfile_)
120 {}
121 
122 
125 (
127 )
128 :
130  flowRate_(ptf.flowRate_.clone()),
131  rhoName_(ptf.rhoName_),
132  rhoInlet_(ptf.rhoInlet_),
133  volumetric_(ptf.volumetric_),
134  extrapolateProfile_(ptf.extrapolateProfile_)
135 {}
136 
137 
140 (
143 )
144 :
146  flowRate_(ptf.flowRate_.clone()),
147  rhoName_(ptf.rhoName_),
148  rhoInlet_(ptf.rhoInlet_),
149  volumetric_(ptf.volumetric_),
150  extrapolateProfile_(ptf.extrapolateProfile_)
151 {}
152 
153 
154 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
155 
156 template<class RhoType>
157 void Foam::flowRateInletVelocityFvPatchVectorField::updateValues
158 (
159  const RhoType& rho
160 )
161 {
162  const scalar t = db().time().timeOutputValue();
163 
164  const vectorField n(patch().nf());
165 
166  if (extrapolateProfile_)
167  {
168  vectorField Up(this->patchInternalField());
169 
170  // Patch normal extrapolated velocity
171  scalarField nUp(n & Up);
172 
173  // Remove the normal component of the extrapolate patch velocity
174  Up -= nUp*n;
175 
176  // Remove any reverse flow
177  nUp = min(nUp, scalar(0));
178 
179  const scalar flowRate = flowRate_->value(t);
180  const scalar estimatedFlowRate = -gSum(rho*(this->patch().magSf()*nUp));
181 
182  if (estimatedFlowRate > 0.5*flowRate)
183  {
184  nUp *= (mag(flowRate)/mag(estimatedFlowRate));
185  }
186  else
187  {
188  nUp -= ((flowRate - estimatedFlowRate)/gSum(rho*patch().magSf()));
189  }
190 
191  // Add the corrected normal component of velocity to the patch velocity
192  Up += nUp*n;
193 
194  // Correct the patch velocity
195  this->operator==(Up);
196  }
197  else
198  {
199  const scalar avgU = -flowRate_->value(t)/gSum(rho*patch().magSf());
200  operator==(avgU*n);
201  }
202 }
203 
204 
206 {
207  if (updated())
208  {
209  return;
210  }
211 
212  if (volumetric_ || rhoName_ == "none")
213  {
214  updateValues(one{});
215  }
216  else
217  {
218  // Mass flow-rate
219  if (db().foundObject<volScalarField>(rhoName_))
220  {
221  const fvPatchField<scalar>& rhop =
222  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
223 
224  updateValues(rhop);
225  }
226  else
227  {
228  // Use constant density
229  if (rhoInlet_ < 0)
230  {
232  << "Did not find registered density field " << rhoName_
233  << " and no constant density 'rhoInlet' specified"
234  << exit(FatalError);
235  }
236 
237  updateValues(rhoInlet_);
238  }
239  }
240 
241  fixedValueFvPatchVectorField::updateCoeffs();
242 }
243 
244 
246 {
248  flowRate_->writeData(os);
249  if (!volumetric_)
250  {
251  os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
252  os.writeEntryIfDifferent<scalar>("rhoInlet", -VGREAT, rhoInlet_);
253  }
254  if (extrapolateProfile_)
255  {
256  os.writeEntry("extrapolateProfile", extrapolateProfile_);
257  }
258  writeEntry("value", os);
259 }
260 
261 
262 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
263 
264 namespace Foam
265 {
267  (
270  );
271 }
272 
273 
274 // ************************************************************************* //
flowRateInletVelocityFvPatchVectorField.H
Foam::fvPatchField
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: volSurfaceMapping.H:51
Foam::flowRateInletVelocityFvPatchVectorField::write
virtual void write(Ostream &) const
Write.
Definition: flowRateInletVelocityFvPatchVectorField.C:245
volFields.H
Foam::fvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:384
Foam::Ostream::writeEntryIfDifferent
Ostream & writeEntryIfDifferent(const word &key, const T &value1, const T &value2)
Write a keyword/value entry only when the two values differ.
Definition: Ostream.H:248
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::Switch
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:77
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::one
A class representing the concept of 1 (one) that can be used to avoid manipulating objects known to b...
Definition: one.H:61
Foam::FatalIOError
IOerror FatalIOError
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
rho
rho
Definition: readInitialConditions.H:88
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::fixedValueFvPatchField
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
Definition: fixedValueFvPatchField.H:80
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:54
n
label n
Definition: TABSMDCalcMethod2.H:31
one.H
Foam::Field< vector >
Foam::operator==
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::flowRateInletVelocityFvPatchVectorField::flowRateInletVelocityFvPatchVectorField
flowRateInletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
Definition: flowRateInletVelocityFvPatchVectorField.C:39
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Switch.H
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
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::flowRateInletVelocityFvPatchVectorField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: flowRateInletVelocityFvPatchVectorField.C:205
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::flowRateInletVelocityFvPatchVectorField
Velocity inlet boundary condition either correcting the extrapolated velocity or creating a uniform v...
Definition: flowRateInletVelocityFvPatchVectorField.H:145
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
FatalIOErrorInFunction
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:473
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
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::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::stringOps::evaluate
string evaluate(label fieldWidth, const std::string &s, size_t pos=0, size_t len=std::string::npos)
String evaluation with specified (positive, non-zero) field width.
Definition: stringOpsEvaluate.C:37