totalPressureFvPatchScalarField.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-2019 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 "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
31 
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const fvPatch& p,
39 )
40 :
41  fixedValueFvPatchScalarField(p, iF),
42  UName_("U"),
43  phiName_("phi"),
44  rhoName_("rho"),
45  psiName_("none"),
46  gamma_(0.0),
47  p0_(p.size(), 0.0)
48 {}
49 
50 
52 (
53  const fvPatch& p,
55  const dictionary& dict
56 )
57 :
58  fixedValueFvPatchScalarField(p, iF, dict, false),
59  UName_(dict.lookupOrDefault<word>("U", "U")),
60  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
61  rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
62  psiName_(dict.lookupOrDefault<word>("psi", "none")),
63  gamma_(psiName_ != "none" ? readScalar(dict.lookup("gamma")) : 1),
64  p0_("p0", dict, p.size())
65 {
66  if (dict.found("value"))
67  {
69  (
70  scalarField("value", dict, p.size())
71  );
72  }
73  else
74  {
76  }
77 }
78 
79 
81 (
83  const fvPatch& p,
85  const fvPatchFieldMapper& mapper
86 )
87 :
88  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
89  UName_(ptf.UName_),
90  phiName_(ptf.phiName_),
91  rhoName_(ptf.rhoName_),
92  psiName_(ptf.psiName_),
93  gamma_(ptf.gamma_),
94  p0_(mapper(ptf.p0_))
95 {}
96 
97 
99 (
101 )
102 :
103  fixedValueFvPatchScalarField(tppsf),
104  UName_(tppsf.UName_),
105  phiName_(tppsf.phiName_),
106  rhoName_(tppsf.rhoName_),
107  psiName_(tppsf.psiName_),
108  gamma_(tppsf.gamma_),
109  p0_(tppsf.p0_)
110 {}
111 
112 
114 (
115  const totalPressureFvPatchScalarField& tppsf,
117 )
118 :
119  fixedValueFvPatchScalarField(tppsf, iF),
120  UName_(tppsf.UName_),
121  phiName_(tppsf.phiName_),
122  rhoName_(tppsf.rhoName_),
123  psiName_(tppsf.psiName_),
124  gamma_(tppsf.gamma_),
125  p0_(tppsf.p0_)
126 {}
127 
128 
129 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
130 
132 (
133  const fvPatchFieldMapper& m
134 )
135 {
136  fixedValueFvPatchScalarField::autoMap(m);
137  m(p0_, p0_);
138 }
139 
140 
142 (
143  const fvPatchScalarField& ptf,
144  const labelList& addr
145 )
146 {
147  fixedValueFvPatchScalarField::rmap(ptf, addr);
148 
149  const totalPressureFvPatchScalarField& tiptf =
150  refCast<const totalPressureFvPatchScalarField>(ptf);
151 
152  p0_.rmap(tiptf.p0_, addr);
153 }
154 
155 
157 (
158  const scalarField& p0p,
159  const vectorField& Up
160 )
161 {
162  if (updated())
163  {
164  return;
165  }
166 
167  const fvsPatchField<scalar>& phip =
168  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
169 
170  if (internalField().dimensions() == dimPressure)
171  {
172  if (psiName_ == "none")
173  {
174  // Variable density and low-speed compressible flow
175 
176  const fvPatchField<scalar>& rho =
177  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
178 
179  operator==(p0p - 0.5*rho*(1.0 - pos0(phip))*magSqr(Up));
180  }
181  else
182  {
183  // High-speed compressible flow
184 
185  const fvPatchField<scalar>& psip =
186  patch().lookupPatchField<volScalarField, scalar>(psiName_);
187 
188  if (gamma_ > 1)
189  {
190  scalar gM1ByG = (gamma_ - 1)/gamma_;
191 
192  operator==
193  (
194  p0p
195  /pow
196  (
197  (1.0 + 0.5*psip*gM1ByG*(1.0 - pos0(phip))*magSqr(Up)),
198  1.0/gM1ByG
199  )
200  );
201  }
202  else
203  {
204  operator==(p0p/(1.0 + 0.5*psip*(1.0 - pos0(phip))*magSqr(Up)));
205  }
206  }
207 
208  }
209  else if (internalField().dimensions() == dimPressure/dimDensity)
210  {
211  // Incompressible flow
212  operator==(p0p - 0.5*(1.0 - pos0(phip))*magSqr(Up));
213  }
214  else
215  {
217  << " Incorrect pressure dimensions " << internalField().dimensions()
218  << nl
219  << " Should be " << dimPressure
220  << " for compressible/variable density flow" << nl
221  << " or " << dimPressure/dimDensity
222  << " for incompressible flow," << nl
223  << " on patch " << this->patch().name()
224  << " of field " << this->internalField().name()
225  << " in file " << this->internalField().objectPath()
226  << exit(FatalError);
227  }
228 
229  fixedValueFvPatchScalarField::updateCoeffs();
230 }
231 
232 
234 {
236  (
237  p0(),
238  patch().lookupPatchField<volVectorField, vector>(UName())
239  );
240 }
241 
242 
244 {
246  writeEntryIfDifferent<word>(os, "U", "U", UName_);
247  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
248  writeEntry(os, "rho", rhoName_);
249  writeEntry(os, "psi", psiName_);
250  writeEntry(os, "gamma", gamma_);
251  writeEntry(os, "p0", p0_);
252  writeEntry(os, "value", *this);
253 }
254 
255 
256 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
257 
258 namespace Foam
259 {
261  (
264  );
265 }
266 
267 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
Foam::surfaceFields.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:438
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
const word & UName() const
Return the name of the velocity field.
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:280
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
const dimensionSet dimPressure
virtual label size() const
Return size.
Definition: fvPatch.H:155
dimensioned< scalar > magSqr(const dimensioned< Type > &)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
volScalarField scalarField(fieldObject, mesh)
dimensionedScalar pos0(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:265
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
const dimensionSet dimDensity
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
totalPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const scalarField & p0() const
Return the total pressure.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:583
This boundary condition provides a total pressure condition. Four variants are possible: ...