uniformTotalPressureFvPatchScalarField.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 | Copyright (C) 2011-2016 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
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_()
48 {}
49 
50 
53 (
54  const fvPatch& p,
56  const dictionary& dict
57 )
58 :
59  fixedValueFvPatchScalarField(p, iF),
60  UName_(dict.lookupOrDefault<word>("U", "U")),
61  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
62  rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
63  psiName_(dict.lookupOrDefault<word>("psi", "none")),
64  gamma_(psiName_ != "none" ? readScalar(dict.lookup("gamma")) : 1),
65  p0_(Function1<scalar>::New("p0", dict))
66 {
67  if (dict.found("value"))
68  {
70  (
71  scalarField("value", dict, p.size())
72  );
73  }
74  else
75  {
76  const scalar t = this->db().time().timeOutputValue();
77  fvPatchScalarField::operator==(p0_->value(t));
78  }
79 }
80 
81 
84 (
86  const fvPatch& p,
88  const fvPatchFieldMapper& mapper
89 )
90 :
91  fixedValueFvPatchScalarField(p, iF), // Don't map
92  UName_(ptf.UName_),
93  phiName_(ptf.phiName_),
94  rhoName_(ptf.rhoName_),
95  psiName_(ptf.psiName_),
96  gamma_(ptf.gamma_),
97  p0_(ptf.p0_, false)
98 {
99  patchType() = ptf.patchType();
100 
101  // Set the patch pressure to the current total pressure
102  // This is not ideal but avoids problems with the creation of patch faces
103  const scalar t = this->db().time().timeOutputValue();
104  fvPatchScalarField::operator==(p0_->value(t));
105 }
106 
107 
110 (
112 )
113 :
114  fixedValueFvPatchScalarField(ptf),
115  UName_(ptf.UName_),
116  phiName_(ptf.phiName_),
117  rhoName_(ptf.rhoName_),
118  psiName_(ptf.psiName_),
119  gamma_(ptf.gamma_),
120  p0_(ptf.p0_, false)
121 {}
122 
123 
126 (
129 )
130 :
131  fixedValueFvPatchScalarField(ptf, iF),
132  UName_(ptf.UName_),
133  phiName_(ptf.phiName_),
134  rhoName_(ptf.rhoName_),
135  psiName_(ptf.psiName_),
136  gamma_(ptf.gamma_),
137  p0_(ptf.p0_, false)
138 {}
139 
140 
141 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
142 
144 (
145  const vectorField& Up
146 )
147 {
148  if (updated())
149  {
150  return;
151  }
152 
153  scalar p0 = p0_->value(this->db().time().timeOutputValue());
154 
155  const fvsPatchField<scalar>& phip =
156  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
157 
158  if (internalField().dimensions() == dimPressure)
159  {
160  if (psiName_ == "none")
161  {
162  // Variable density and low-speed compressible flow
163 
164  const fvPatchField<scalar>& rho =
165  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
166 
167  operator==(p0 - 0.5*rho*(1.0 - pos(phip))*magSqr(Up));
168  }
169  else
170  {
171  // High-speed compressible flow
172 
173  const fvPatchField<scalar>& psip =
174  patch().lookupPatchField<volScalarField, scalar>(psiName_);
175 
176  if (gamma_ > 1)
177  {
178  scalar gM1ByG = (gamma_ - 1)/gamma_;
179 
180  operator==
181  (
182  p0
183  /pow
184  (
185  (1.0 + 0.5*psip*gM1ByG*(1.0 - pos(phip))*magSqr(Up)),
186  1.0/gM1ByG
187  )
188  );
189  }
190  else
191  {
192  operator==(p0/(1.0 + 0.5*psip*(1.0 - pos(phip))*magSqr(Up)));
193  }
194  }
195 
196  }
197  else if (internalField().dimensions() == dimPressure/dimDensity)
198  {
199  // Incompressible flow
200  operator==(p0 - 0.5*(1.0 - pos(phip))*magSqr(Up));
201  }
202  else
203  {
205  << " Incorrect pressure dimensions " << internalField().dimensions()
206  << nl
207  << " Should be " << dimPressure
208  << " for compressible/variable density flow" << nl
209  << " or " << dimPressure/dimDensity
210  << " for incompressible flow," << nl
211  << " on patch " << this->patch().name()
212  << " of field " << this->internalField().name()
213  << " in file " << this->internalField().objectPath()
214  << exit(FatalError);
215  }
216 
217  fixedValueFvPatchScalarField::updateCoeffs();
218 }
219 
220 
222 {
223  updateCoeffs(patch().lookupPatchField<volVectorField, vector>(UName_));
224 }
225 
226 
228 {
230  writeEntryIfDifferent<word>(os, "U", "U", UName_);
231  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
232  os.writeKeyword("rho") << rhoName_ << token::END_STATEMENT << nl;
233  os.writeKeyword("psi") << psiName_ << token::END_STATEMENT << nl;
234  os.writeKeyword("gamma") << gamma_ << token::END_STATEMENT << nl;
235  p0_->writeData(os);
236  writeEntry("value", os);
237 }
238 
239 
240 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
241 
242 namespace Foam
243 {
245  (
248  );
249 }
250 
251 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
Foam::surfaceFields.
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: Function1.H:53
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
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
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:65
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
dimensionedScalar pos(const dimensionedScalar &ds)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
virtual label size() const
Return size.
Definition: fvPatch.H:161
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
const dimensionSet dimPressure
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)
uniformTotalPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
static const char nl
Definition: Ostream.H:262
This boundary condition provides a time-varying form of the uniform total pressure boundary condition...
const dimensionSet dimDensity
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
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.
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:363
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:451