matchedFlowRateOutletVelocityFvPatchVectorField.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-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 
27 #include "volFields.H"
28 #include "one.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
35 (
36  const fvPatch& p,
38 )
39 :
41  inletPatchName_(),
42  volumetric_(false),
43  rhoName_("rho")
44 {}
45 
46 
49 (
50  const fvPatch& p,
52  const dictionary& dict
53 )
54 :
56  inletPatchName_(dict.lookup("inletPatch")),
57  volumetric_(dict.lookupOrDefault("volumetric", true))
58 {
59  if (volumetric_)
60  {
61  rhoName_ = "none";
62  }
63  else
64  {
65  rhoName_ = word(dict.lookupOrDefault<word>("rho", "rho"));
66  }
67 
68  // Value field require if mass based
69  if (dict.found("value"))
70  {
72  (
73  vectorField("value", dict, p.size())
74  );
75  }
76  else
77  {
78  evaluate(Pstream::commsTypes::blocking);
79  }
80 }
81 
82 
85 (
87  const fvPatch& p,
89  const fvPatchFieldMapper& mapper
90 )
91 :
92  fixedValueFvPatchField<vector>(ptf, p, iF, mapper),
93  inletPatchName_(ptf.inletPatchName_),
94  volumetric_(ptf.volumetric_),
95  rhoName_(ptf.rhoName_)
96 {}
97 
98 
101 (
104 )
105 :
107  inletPatchName_(ptf.inletPatchName_),
108  volumetric_(ptf.volumetric_),
109  rhoName_(ptf.rhoName_)
110 {}
111 
112 
113 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
114 
115 template<class RhoType>
116 void Foam::matchedFlowRateOutletVelocityFvPatchVectorField::updateValues
117 (
118  const label inletPatchID,
119  const RhoType& rhoOutlet,
120  const RhoType& rhoInlet
121 )
122 {
123  const fvPatch& p = patch();
124  const fvPatch& inletPatch = p.boundaryMesh()[inletPatchID];
125 
126  const vectorField n(p.nf());
127 
128  // Extrapolate patch velocity
129  vectorField Up(patchInternalField());
130 
131  // Patch normal extrapolated velocity
132  scalarField nUp(n & Up);
133 
134  // Remove the normal component of the extrapolate patch velocity
135  Up -= nUp*n;
136 
137  // Remove any reverse flow
138  nUp = max(nUp, scalar(0));
139 
140  // Lookup non-const access to velocity field
142  (
143  const_cast<volVectorField&>
144  (
145  dynamic_cast<const volVectorField&>(internalField())
146  )
147  );
148 
149  // Get the corresponding inlet velocity patch field
150  fvPatchVectorField& inletPatchU = U.boundaryFieldRef()[inletPatchID];
151 
152  // Ensure that the corresponding inlet velocity patch field is up-to-date
153  inletPatchU.updateCoeffs();
154 
155  // Calculate the inlet patch flow rate
156  const scalar flowRate = -gSum(rhoInlet*(inletPatch.Sf() & inletPatchU));
157 
158  // Calculate the extrapolated outlet patch flow rate
159  const scalar estimatedFlowRate = gSum(rhoOutlet*(patch().magSf()*nUp));
160 
161  if (estimatedFlowRate/flowRate > 0.5)
162  {
163  nUp *= (mag(flowRate)/mag(estimatedFlowRate));
164  }
165  else
166  {
167  nUp += ((flowRate - estimatedFlowRate)/gSum(rhoOutlet*patch().magSf()));
168  }
169 
170  // Add the corrected normal component of velocity to the patch velocity
171  Up += nUp*n;
172 
173  // Correct the patch velocity
174  operator==(Up);
175 }
176 
177 
179 {
180  if (updated())
181  {
182  return;
183  }
184 
185  // Find corresponding inlet patch
186  const label inletPatchID =
187  patch().patch().boundaryMesh().findPatchID(inletPatchName_);
188 
189  if (inletPatchID < 0)
190  {
192  << "Unable to find inlet patch " << inletPatchName_
193  << exit(FatalError);
194  }
195 
196  if (volumetric_)
197  {
198  updateValues(inletPatchID, one(), one());
199  }
200  else
201  {
202  // Mass flow-rate
203  if (db().foundObject<volScalarField>(rhoName_))
204  {
205  const volScalarField& rho = db().lookupObject<volScalarField>
206  (
207  rhoName_
208  );
209 
210  updateValues
211  (
212  inletPatchID,
213  rho.boundaryField()[patch().index()],
214  rho.boundaryField()[inletPatchID]
215  );
216  }
217  else
218  {
220  << "Cannot find density field " << rhoName_ << exit(FatalError);
221  }
222  }
223 
224  fixedValueFvPatchVectorField::updateCoeffs();
225 }
226 
227 
229 (
230  Ostream& os
231 ) const
232 {
234  writeEntry(os, "inletPatch", inletPatchName_);
235  if (!volumetric_)
236  {
237  writeEntry(os, "volumetric", volumetric_);
238  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
239  }
240  writeEntry(os, "value", *this);
241 }
242 
243 
244 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
245 
246 namespace Foam
247 {
249  (
252  );
253 }
254 
255 
256 // ************************************************************************* //
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:180
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void evaluate(GeometricField< Type, PatchField, GeoMesh > &result, const Function1< Type > &func, const GeometricField< Type, PatchField, GeoMesh > &x)
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const Boundary & boundaryField() const
Return const-reference to the boundary field.
tmp< vectorField > nf() const
Return face normals.
Definition: fvPatch.C:130
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
matchedFlowRateOutletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
volVectorField vectorField(fieldObject, mesh)
Velocity outlet boundary condition which corrects the extrapolated velocity to match the flow rate of...
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:260
Macros for easy insertion into run-time selection tables.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Type gSum(const FieldField< Field, Type > &f)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
const vectorField & Sf() const
Return face area vectors.
Definition: fvPatch.C:136
U
Definition: pEqn.H:72
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:209
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
label n
volScalarField & p
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
A class representing the concept of 1 (scalar(1)) used to avoid unnecessary manipulations for objects...
Definition: one.H:50