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-2018 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 (
103 )
104 :
106  inletPatchName_(ptf.inletPatchName_),
107  volumetric_(ptf.volumetric_),
108  rhoName_(ptf.rhoName_)
109 {}
110 
111 
114 (
117 )
118 :
120  inletPatchName_(ptf.inletPatchName_),
121  volumetric_(ptf.volumetric_),
122  rhoName_(ptf.rhoName_)
123 {}
124 
125 
126 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
127 
128 template<class RhoType>
129 void Foam::matchedFlowRateOutletVelocityFvPatchVectorField::updateValues
130 (
131  const label inletPatchID,
132  const RhoType& rhoOutlet,
133  const RhoType& rhoInlet
134 )
135 {
136  const fvPatch& p = patch();
137  const fvPatch& inletPatch = p.boundaryMesh()[inletPatchID];
138 
139  const vectorField n(p.nf());
140 
141  // Extrapolate patch velocity
142  vectorField Up(patchInternalField());
143 
144  // Patch normal extrapolated velocity
145  scalarField nUp(n & Up);
146 
147  // Remove the normal component of the extrapolate patch velocity
148  Up -= nUp*n;
149 
150  // Remove any reverse flow
151  nUp = max(nUp, scalar(0));
152 
153  // Lookup non-const access to velocity field
155  (
156  const_cast<volVectorField&>
157  (
158  dynamic_cast<const volVectorField&>(internalField())
159  )
160  );
161 
162  // Get the corresponding inlet velocity patch field
163  fvPatchVectorField& inletPatchU = U.boundaryFieldRef()[inletPatchID];
164 
165  // Ensure that the corresponding inlet velocity patch field is up-to-date
166  inletPatchU.updateCoeffs();
167 
168  // Calculate the inlet patch flow rate
169  const scalar flowRate = -gSum(rhoInlet*(inletPatch.Sf() & inletPatchU));
170 
171  // Calculate the extrapolated outlet patch flow rate
172  const scalar estimatedFlowRate = gSum(rhoOutlet*(patch().magSf()*nUp));
173 
174  if (estimatedFlowRate/flowRate > 0.5)
175  {
176  nUp *= (mag(flowRate)/mag(estimatedFlowRate));
177  }
178  else
179  {
180  nUp += ((flowRate - estimatedFlowRate)/gSum(rhoOutlet*patch().magSf()));
181  }
182 
183  // Add the corrected normal component of velocity to the patch velocity
184  Up += nUp*n;
185 
186  // Correct the patch velocity
187  operator==(Up);
188 }
189 
190 
192 {
193  if (updated())
194  {
195  return;
196  }
197 
198  // Find corresponding inlet patch
199  const label inletPatchID =
200  patch().patch().boundaryMesh().findPatchID(inletPatchName_);
201 
202  if (inletPatchID < 0)
203  {
205  << "Unable to find inlet patch " << inletPatchName_
206  << exit(FatalError);
207  }
208 
209  if (volumetric_)
210  {
211  updateValues(inletPatchID, one(), one());
212  }
213  else
214  {
215  // Mass flow-rate
216  if (db().foundObject<volScalarField>(rhoName_))
217  {
218  const volScalarField& rho = db().lookupObject<volScalarField>
219  (
220  rhoName_
221  );
222 
223  updateValues
224  (
225  inletPatchID,
226  rho.boundaryField()[patch().index()],
227  rho.boundaryField()[inletPatchID]
228  );
229  }
230  else
231  {
233  << "Cannot find density field " << rhoName_ << exit(FatalError);
234  }
235  }
236 
237  fixedValueFvPatchVectorField::updateCoeffs();
238 }
239 
240 
242 (
243  Ostream& os
244 ) const
245 {
247  os.writeKeyword("inletPatch")
248  << inletPatchName_ << token::END_STATEMENT << nl;
249  if (!volumetric_)
250  {
251  os.writeKeyword("volumetric")
252  << volumetric_ << token::END_STATEMENT << nl;
253  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
254  }
255  writeEntry("value", os);
256 }
257 
258 
259 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
260 
261 namespace Foam
262 {
264  (
267  );
268 }
269 
270 
271 // ************************************************************************* //
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:185
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
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const Boundary & boundaryField() const
Return const-reference to the boundary field.
tmp< vectorField > nf() const
Return face normals.
Definition: fvPatch.C:124
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
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:362
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:53
static const char nl
Definition: Ostream.H:265
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
const vectorField & Sf() const
Return face area vectors.
Definition: fvPatch.C:130
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:311
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