pressureDirectedInletOutletVelocityFvPatchVectorField.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-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"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
36 (
37  const fvPatch& p,
39 )
40 :
41  mixedFvPatchVectorField(p, iF),
42  phiName_("phi"),
43  rhoName_("rho"),
44  inletDir_(p.size())
45 {
46  refValue() = *this;
47  refGrad() = Zero;
48  valueFraction() = 0.0;
49 }
50 
51 
54 (
56  const fvPatch& p,
58  const fvPatchFieldMapper& mapper
59 )
60 :
61  mixedFvPatchVectorField(ptf, p, iF, mapper),
62  phiName_(ptf.phiName_),
63  rhoName_(ptf.rhoName_),
64  inletDir_(ptf.inletDir_, mapper)
65 {}
66 
67 
70 (
71  const fvPatch& p,
73  const dictionary& dict
74 )
75 :
76  mixedFvPatchVectorField(p, iF),
77  phiName_(dict.lookupOrDefault<word>("phi", "phi")),
78  rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
79  inletDir_("inletDirection", dict, p.size())
80 {
81  fvPatchVectorField::operator=(vectorField("value", dict, p.size()));
82  refValue() = *this;
83  refGrad() = Zero;
84  valueFraction() = 0.0;
85 }
86 
87 
90 (
92 )
93 :
94  mixedFvPatchVectorField(pivpvf),
95  phiName_(pivpvf.phiName_),
96  rhoName_(pivpvf.rhoName_),
97  inletDir_(pivpvf.inletDir_)
98 {}
99 
100 
103 (
106 )
107 :
108  mixedFvPatchVectorField(pivpvf, iF),
109  phiName_(pivpvf.phiName_),
110  rhoName_(pivpvf.rhoName_),
111  inletDir_(pivpvf.inletDir_)
112 {}
113 
114 
115 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
116 
118 (
119  const fvPatchFieldMapper& m
120 )
121 {
122  mixedFvPatchVectorField::autoMap(m);
123  inletDir_.autoMap(m);
124 }
125 
126 
128 (
129  const fvPatchVectorField& ptf,
130  const labelList& addr
131 )
132 {
133  mixedFvPatchVectorField::rmap(ptf, addr);
134 
136  refCast<const pressureDirectedInletOutletVelocityFvPatchVectorField>
137  (ptf);
138 
139  inletDir_.rmap(tiptf.inletDir_, addr);
140 }
141 
142 
144 {
145  if (updated())
146  {
147  return;
148  }
149 
150  const surfaceScalarField& phi =
151  db().lookupObject<surfaceScalarField>(phiName_);
152 
153  const fvsPatchField<scalar>& phip =
154  patch().patchField<surfaceScalarField, scalar>(phi);
155 
156  tmp<vectorField> n = patch().nf();
157  tmp<scalarField> ndmagS = (n & inletDir_)*patch().magSf();
158 
159  if (phi.dimensions() == dimVelocity*dimArea)
160  {
161  refValue() = inletDir_*phip/ndmagS;
162  }
163  else if (phi.dimensions() == dimDensity*dimVelocity*dimArea)
164  {
165  const fvPatchField<scalar>& rhop =
166  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
167 
168  refValue() = inletDir_*phip/(rhop*ndmagS);
169  }
170  else
171  {
173  << "dimensions of phi are not correct"
174  << "\n on patch " << this->patch().name()
175  << " of field " << this->internalField().name()
176  << " in file " << this->internalField().objectPath()
177  << exit(FatalError);
178  }
179 
180  valueFraction() = 1.0 - pos0(phip);
181 
182  mixedFvPatchVectorField::updateCoeffs();
183 }
184 
185 
187 (
188  Ostream& os
189 ) const
190 {
192  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
193  writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_);
194  inletDir_.writeEntry("inletDirection", os);
195  writeEntry("value", os);
196 }
197 
198 
199 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
200 
201 void Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::operator=
202 (
203  const fvPatchField<vector>& pvf
204 )
205 {
207  (
208  valueFraction()*(inletDir_*(inletDir_ & pvf))
209  + (1 - valueFraction())*pvf
210  );
211 }
212 
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 namespace Foam
217 {
219  (
222  );
223 }
224 
225 // ************************************************************************* //
Foam::surfaceFields.
surfaceScalarField & phi
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
pressureDirectedInletOutletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
#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:66
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:726
volVectorField vectorField(fieldObject, mesh)
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Macros for easy insertion into run-time selection tables.
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
const dimensionSet & dimensions() const
Return dimensions.
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:91
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual label size() const
Return size.
Definition: fvPatch.H:161
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimDensity
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
This velocity inlet/outlet boundary condition is applied to pressure boundaries where the pressure is...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
label n
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:65
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:57
Namespace for OpenFOAM.
const dimensionSet dimVelocity