outletMachNumberPressureFvPatchScalarField.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) 2018-2020 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
30#include "fvPatchFieldMapper.H"
31#include "volFields.H"
32#include "surfaceFields.H"
33#include "fluidThermo.H"
34#include "constants.H"
35
36// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37
38Foam::outletMachNumberPressureFvPatchScalarField::
39outletMachNumberPressureFvPatchScalarField
40(
41 const fvPatch& p,
43)
44:
45 fixedValueFvPatchScalarField(p, iF),
46 M_(1),
47 pBack_(0.0),
48 c1_(0.0),
49 A1_(0.0),
50 phiName_("phi"),
51 rhoName_("rho"),
52 UName_("U"),
53 choked_(false),
54 relax_(0.0)
55{}
56
57
58Foam::outletMachNumberPressureFvPatchScalarField::
59outletMachNumberPressureFvPatchScalarField
60(
61 const fvPatch& p,
63 const dictionary& dict
64)
65:
66 fixedValueFvPatchScalarField(p, iF, dict),
67 M_(dict.getOrDefault<scalar>("M", 0)),
68 pBack_(dict.get<scalar>("pBack")),
69 c1_(dict.getOrDefault<scalar>("c1", 0)),
70 A1_(dict.getOrDefault<scalar>("A1", 0)),
71 phiName_(dict.getOrDefault<word>("phi", "phi")),
72 rhoName_(dict.getOrDefault<word>("rho", "rho")),
73 UName_(dict.getOrDefault<word>("U", "U")),
74 choked_(dict.get<Switch>("choked")),
75 relax_(dict.getOrDefault<scalar>("relax", 0))
76{}
77
78
79Foam::outletMachNumberPressureFvPatchScalarField::
80outletMachNumberPressureFvPatchScalarField
81(
83 const fvPatch& p,
85 const fvPatchFieldMapper& mapper
86)
87:
88 fixedValueFvPatchScalarField(ptf, p, iF, mapper),
89 M_(ptf.M_),
90 pBack_(ptf.pBack_),
91 c1_(ptf.c1_),
92 A1_(ptf.A1_),
93 phiName_(ptf.phiName_),
94 rhoName_(ptf.rhoName_),
95 UName_(ptf.UName_),
96 choked_(ptf.choked_),
97 relax_(ptf.relax_)
98{}
99
100
101Foam::outletMachNumberPressureFvPatchScalarField::
102outletMachNumberPressureFvPatchScalarField
103(
105)
106:
107 fixedValueFvPatchScalarField(tppsf),
108 M_(tppsf.M_),
109 pBack_(tppsf.pBack_),
110 c1_(tppsf.c1_),
111 A1_(tppsf.A1_),
112 phiName_(tppsf.phiName_),
113 rhoName_(tppsf.rhoName_),
114 UName_(tppsf.UName_),
115 choked_(tppsf.choked_),
116 relax_(tppsf.relax_)
117{}
118
119
120Foam::outletMachNumberPressureFvPatchScalarField::
121outletMachNumberPressureFvPatchScalarField
122(
125)
126:
127 fixedValueFvPatchScalarField(tppsf, iF),
128 M_(tppsf.M_),
129 pBack_(tppsf.pBack_),
130 c1_(tppsf.c1_),
131 A1_(tppsf.A1_),
132 phiName_(tppsf.phiName_),
133 rhoName_(tppsf.rhoName_),
134 UName_(tppsf.UName_),
135 choked_(tppsf.choked_),
136 relax_(tppsf.relax_)
137{}
138
139
140// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
141
143{
144 if (updated())
145 {
146 return;
147 }
148
149 const volScalarField& p =
150 db().lookupObject<volScalarField>
151 (
152 this->internalField().name()
153 );
154
155 const label patchi = patch().index();
156
157 const scalarField pb(p.oldTime().boundaryField()[patchi]);
158
160 patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
161
162 // Calculate the current mass flow rate
164 {
166 <<"phi is not a mass flux." << exit(FatalError);
167 }
168
169 const fluidThermo* thermoPtr =
170 db().findObject<fluidThermo>(basicThermo::dictName);
171
172 const volVectorField& U = db().lookupObject<volVectorField>(UName_);
173
174 vectorField Ub(U.boundaryField()[patchi]);
175 const vectorField UbOld(U.oldTime().boundaryField()[patchi]);
176
177 // relax U
178 Ub = relax_*UbOld + (1 - relax_)*Ub;
179
180 const scalarField gamma(thermoPtr->gamma()().boundaryField()[patchi]);
181
183 patch().lookupPatchField<volScalarField, scalar>(rhoName_);
184
185 const scalarField Mb(mag(Ub)/sqrt(gamma*pb/rho));
186
187 const scalarField ptot
188 (
189 pb*(pow(1 + (gamma-1)/2*sqr(gAverage(Mb)), gamma/(gamma-1)))
190 );
191
192 scalarField M(patch().size(), 1.0);
193
194 if (choked_)
195 {
196 if (M_ > 0.0)
197 {
198 M = M_;
199 }
200 else
201 {
202 FatalErrorInFunction <<" Mach number is lower than zero" << endl
203 << "Pelase specify M in the dictionary"
204 << exit(FatalError);
205 }
206 }
207 else
208 {
209
210 if (A1_ == 0.0 || c1_ == 0.0)
211 {
212 FatalErrorInFunction <<" Please enter values for A1 and c1" << endl
213 << exit(FatalError);
214 }
215
216 const scalarField r(pBack_/ptot);
217 const scalar area = gSum(mag(patch().Sf()));
218 M =
219 A1_/(c1_*area)
220 *sqrt(2/(gamma-1)*(pow(r, 2/gamma) - pow(r, (gamma+1)/gamma)));
221
222 forAll(M, i)
223 {
224 if (M[i] < 0 || r[i] >=1)
225 {
226 WarningInFunction <<" Mach number is lower than zero" << endl
227 << "or pBack/ptot ratio is larger then one"
228 << endl;
229 }
230 }
231 }
232
233 scalarField pbNew
234 (
235 ptot/(pow(1 + (gamma-1)/2*sqr(gAverage(M)), gamma/(gamma-1)))
236 );
237
238 // relax pressure
239 pbNew = relax_*pb + (1 -relax_)*pbNew;
240
241 operator==(pbNew);
242
243 fixedValueFvPatchScalarField::updateCoeffs();
244}
245
246
248{
250 os.writeEntry("pBack", pBack_);
251 os.writeEntryIfDifferent<scalar>("c1", 0, c1_);
252 os.writeEntryIfDifferent<scalar>("A1", 0, A1_);
253 os.writeEntry("choked", choked_);
254 os.writeEntryIfDifferent<scalar>("relax", 0, relax_);
255
256 os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
257 os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
258 os.writeEntryIfDifferent<word>("U", "U", UName_);
259 os.writeEntryIfDifferent<scalar>("M", 0, M_);
260
261 writeEntry("value", os);
262}
263
264
265// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
266
267namespace Foam
268{
270 (
273 );
274}
275
276// ************************************************************************* //
#define M(I)
Macros for easy insertion into run-time selection tables.
surfaceScalarField & phi
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return dimensions.
const GeometricField< Type, PatchField, GeoMesh > & oldTime() const
Return old time field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
const Internal & internalField() const
Return a const-reference to the dimensioned internal field.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
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:251
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
static const word dictName
Definition: basicThermo.H:256
virtual tmp< volScalarField > gamma() const =0
Gamma = Cp/Cv [].
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Fundamental fluid thermodynamic properties.
Definition: fluidThermo.H:56
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:79
This boundary condition maintains a certain subsonic Mach number at an outlet patch by dynamically ad...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
volScalarField & p
const scalar gamma
Definition: EEqn.H:9
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
#define makePatchTypeField(PatchTypeField, typePatchTypeField)
Definition: fvPatchField.H:676
#define WarningInFunction
Report a warning using Foam::Warning.
Namespace for OpenFOAM.
Type gSum(const FieldField< Field, Type > &f)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
const dimensionSet dimVelocity
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Type gAverage(const FieldField< Field, Type > &f)
error FatalError
const dimensionSet dimDensity
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
Foam::surfaceFields.