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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2017-2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31#include "fvPatchFieldMapper.H"
32#include "volFields.H"
33#include "surfaceFields.H"
34
35// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36
37Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::
38pressureDirectedInletOutletVelocityFvPatchVectorField
39(
40 const fvPatch& p,
42)
43:
44 mixedFvPatchVectorField(p, iF),
45 phiName_("phi"),
46 rhoName_("rho"),
47 inletDir_(p.size())
48{
49 refValue() = *this;
50 refGrad() = Zero;
51 valueFraction() = 0.0;
52}
53
54
55Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::
56pressureDirectedInletOutletVelocityFvPatchVectorField
57(
59 const fvPatch& p,
61 const fvPatchFieldMapper& mapper
62)
63:
64 mixedFvPatchVectorField(ptf, p, iF, mapper),
65 phiName_(ptf.phiName_),
66 rhoName_(ptf.rhoName_),
67 inletDir_(ptf.inletDir_, mapper)
68{}
69
70
71Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::
72pressureDirectedInletOutletVelocityFvPatchVectorField
73(
74 const fvPatch& p,
76 const dictionary& dict
77)
78:
79 mixedFvPatchVectorField(p, iF),
80 phiName_(dict.getOrDefault<word>("phi", "phi")),
81 rhoName_(dict.getOrDefault<word>("rho", "rho")),
82 inletDir_("inletDirection", dict, p.size())
83{
84 patchType() = dict.getOrDefault<word>("patchType", word::null);
86 refValue() = *this;
87 refGrad() = Zero;
88 valueFraction() = 0.0;
89}
90
91
92Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::
93pressureDirectedInletOutletVelocityFvPatchVectorField
94(
96)
97:
98 mixedFvPatchVectorField(pivpvf),
99 phiName_(pivpvf.phiName_),
100 rhoName_(pivpvf.rhoName_),
101 inletDir_(pivpvf.inletDir_)
102{}
103
104
105Foam::pressureDirectedInletOutletVelocityFvPatchVectorField::
106pressureDirectedInletOutletVelocityFvPatchVectorField
107(
110)
111:
112 mixedFvPatchVectorField(pivpvf, iF),
113 phiName_(pivpvf.phiName_),
114 rhoName_(pivpvf.rhoName_),
115 inletDir_(pivpvf.inletDir_)
116{}
117
118
119// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
120
122(
123 const fvPatchFieldMapper& m
124)
125{
126 mixedFvPatchVectorField::autoMap(m);
127 inletDir_.autoMap(m);
128}
129
130
132(
133 const fvPatchVectorField& ptf,
134 const labelList& addr
135)
136{
137 mixedFvPatchVectorField::rmap(ptf, addr);
138
140 refCast<const pressureDirectedInletOutletVelocityFvPatchVectorField>
141 (ptf);
142
143 inletDir_.rmap(tiptf.inletDir_, addr);
144}
145
146
148{
149 if (updated())
150 {
151 return;
152 }
153
154 const surfaceScalarField& phi =
155 db().lookupObject<surfaceScalarField>(phiName_);
156
157 const fvsPatchField<scalar>& phip =
158 patch().patchField<surfaceScalarField, scalar>(phi);
159
160 tmp<vectorField> n = patch().nf();
161 tmp<scalarField> ndmagS = (n & inletDir_)*patch().magSf();
162
164 {
165 refValue() = inletDir_*phip/ndmagS;
166 }
168 {
169 const fvPatchField<scalar>& rhop =
170 patch().lookupPatchField<volScalarField, scalar>(rhoName_);
171
172 refValue() = inletDir_*phip/(rhop*ndmagS);
173 }
174 else
175 {
177 << "dimensions of phi are not correct"
178 << "\n on patch " << this->patch().name()
179 << " of field " << this->internalField().name()
180 << " in file " << this->internalField().objectPath()
181 << exit(FatalError);
182 }
183
184 valueFraction() = 1.0 - pos0(phip);
185
186 mixedFvPatchVectorField::updateCoeffs();
187}
188
189
191(
192 Ostream& os
193) const
194{
196 os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
197 os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
198 inletDir_.writeEntry("inletDirection", os);
199 writeEntry("value", os);
200}
201
202
203// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
204
206(
207 const fvPatchField<vector>& pvf
208)
209{
211 (
212 valueFraction()*(inletDir_*(inletDir_ & pvf))
213 + (1 - valueFraction())*pvf
214 );
215}
216
217
218// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
219
220namespace Foam
221{
223 (
226 );
227}
228
229// ************************************************************************* //
label n
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.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
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
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
friend Ostream & operator(Ostream &, const faMatrix< Type > &)
virtual bool write()
Write the output fields.
A FieldMapper for finite-volume patch fields.
virtual void operator=(const UList< vector > &)
Definition: fvPatchField.C:408
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 velocity inlet/outlet boundary condition is applied to velocity boundaries where the pressure is...
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void rmap(const fvPatchVectorField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
static const word null
An empty word.
Definition: word.H:80
volScalarField & p
#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
Namespace for OpenFOAM.
dimensionedScalar pos0(const dimensionedScalar &ds)
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
Field< vector > vectorField
Specialisation of Field<T> for vector.
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
const dimensionSet dimDensity
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dictionary dict
Foam::surfaceFields.