prghTotalPressureFvPatchScalarField.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) 2015-2017 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"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
37 (
38  const fvPatch& p,
40 )
41 :
42  fixedValueFvPatchScalarField(p, iF),
43  UName_("U"),
44  phiName_("phi"),
45  rhoName_("rho"),
46  p0_(p.size(), 0.0)
47 {}
48 
49 
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  p0_("p0", dict, p.size())
63 {
64  if (dict.found("value"))
65  {
66  fvPatchScalarField::operator=
67  (
68  scalarField("value", dict, p.size())
69  );
70  }
71  else
72  {
74  }
75 }
76 
77 
80 (
82  const fvPatch& p,
84  const fvPatchFieldMapper& mapper
85 )
86 :
87  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
88  UName_(ptf.UName_),
89  phiName_(ptf.phiName_),
90  rhoName_(ptf.rhoName_),
91  p0_(ptf.p0_, mapper)
92 {}
93 
94 
97 (
99 )
100 :
101  fixedValueFvPatchScalarField(ptf),
102  UName_(ptf.UName_),
103  phiName_(ptf.phiName_),
104  rhoName_(ptf.rhoName_),
105  p0_(ptf.p0_)
106 {}
107 
108 
111 (
114 )
115 :
116  fixedValueFvPatchScalarField(ptf, iF),
117  UName_(ptf.UName_),
118  phiName_(ptf.phiName_),
119  rhoName_(ptf.rhoName_),
120  p0_(ptf.p0_)
121 {}
122 
123 
124 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
125 
127 (
128  const fvPatchFieldMapper& m
129 )
130 {
131  fixedValueFvPatchScalarField::autoMap(m);
132  p0_.autoMap(m);
133 }
134 
135 
137 (
138  const fvPatchScalarField& ptf,
139  const labelList& addr
140 )
141 {
142  fixedValueFvPatchScalarField::rmap(ptf, addr);
143 
145  refCast<const prghTotalPressureFvPatchScalarField>(ptf);
146 
147  p0_.rmap(tiptf.p0_, addr);
148 }
149 
150 
152 {
153  if (updated())
154  {
155  return;
156  }
157 
158  const scalarField& rhop =
159  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
160 
161  const scalarField& phip =
162  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
163 
164  const vectorField& Up =
165  patch().lookupPatchField<volVectorField, vector>(UName_);
166 
168  db().lookupObject<uniformDimensionedVectorField>("g");
169 
170  const uniformDimensionedScalarField& hRef =
171  db().lookupObject<uniformDimensionedScalarField>("hRef");
172 
173  dimensionedScalar ghRef
174  (
175  mag(g.value()) > SMALL
176  ? g & (cmptMag(g.value())/mag(g.value()))*hRef
177  : dimensionedScalar("ghRef", g.dimensions()*dimLength, 0)
178  );
179 
180  operator==
181  (
182  p0_
183  - 0.5*rhop*(1.0 - pos0(phip))*magSqr(Up)
184  - rhop*((g.value() & patch().Cf()) - ghRef.value())
185  );
186 
187  fixedValueFvPatchScalarField::updateCoeffs();
188 }
189 
190 
192 {
194  writeEntryIfDifferent<word>(os, "U", "U", UName_);
195  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
196  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
197  p0_.writeEntry("p0", os);
198  writeEntry("value", os);
199 }
200 
201 
202 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
203 
204 namespace Foam
205 {
207  (
210  );
211 }
212 
213 // ************************************************************************* //
Foam::surfaceFields.
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:431
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
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
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:726
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
This boundary condition provides static pressure condition for p_rgh, calculated as: ...
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
A class for handling words, derived from string.
Definition: word.H:59
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Foam::fvPatchFieldMapper.
virtual label size() const
Return size.
Definition: fvPatch.H:161
void cmptMag(FieldField< Field, Type > &cf, const FieldField< Field, Type > &f)
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)
prghTotalPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
const dimensionedVector & g
word phiName_
Name of the flux transporting the field.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Namespace for OpenFOAM.