fanPressureFvPatchScalarField.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 "volFields.H"
29 #include "surfaceFields.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35  template<>
36  const char* NamedEnum
37  <
39  2
40  >::names[] =
41  {
42  "in",
43  "out"
44  };
45 }
46 
47 const Foam::NamedEnum
48 <
50  2
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
57 (
58  const fvPatch& p,
60 )
61 :
63  fanCurve_(),
64  direction_(ffdOut)
65 {}
66 
67 
69 (
71  const fvPatch& p,
73  const fvPatchFieldMapper& mapper
74 )
75 :
76  totalPressureFvPatchScalarField(ptf, p, iF, mapper),
77  fanCurve_(ptf.fanCurve_),
78  direction_(ptf.direction_)
79 {}
80 
81 
83 (
84  const fvPatch& p,
86  const dictionary& dict
87 )
88 :
90  fanCurve_(dict),
91  direction_(fanFlowDirectionNames_.read(dict.lookup("direction")))
92 {}
93 
94 
96 (
97  const fanPressureFvPatchScalarField& pfopsf
98 )
99 :
101  fanCurve_(pfopsf.fanCurve_),
102  direction_(pfopsf.direction_)
103 {}
104 
105 
107 (
108  const fanPressureFvPatchScalarField& pfopsf,
110 )
111 :
113  fanCurve_(pfopsf.fanCurve_),
114  direction_(pfopsf.direction_)
115 {}
116 
117 
118 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
119 
121 {
122  if (updated())
123  {
124  return;
125  }
126 
127  // Retrieve flux field
128  const surfaceScalarField& phi =
129  db().lookupObject<surfaceScalarField>(phiName());
130 
131  const fvsPatchField<scalar>& phip =
132  patch().patchField<surfaceScalarField, scalar>(phi);
133 
134  int dir = 2*direction_ - 1;
135 
136  // Average volumetric flow rate
137  scalar volFlowRate = 0;
138 
139  if (phi.dimensions() == dimVelocity*dimArea)
140  {
141  volFlowRate = dir*gSum(phip);
142  }
143  else if (phi.dimensions() == dimVelocity*dimArea*dimDensity)
144  {
145  const scalarField& rhop =
146  patch().lookupPatchField<volScalarField, scalar>(rhoName());
147  volFlowRate = dir*gSum(phip/rhop);
148  }
149  else
150  {
152  << "dimensions of phi are not correct"
153  << "\n on patch " << patch().name()
154  << " of field " << internalField().name()
155  << " in file " << internalField().objectPath() << nl
156  << exit(FatalError);
157  }
158 
159  // Pressure drop for this flow rate
160  const scalar pdFan = fanCurve_(max(volFlowRate, 0.0));
161 
163  (
164  p0() - dir*pdFan,
165  patch().lookupPatchField<volVectorField, vector>(UName())
166  );
167 }
168 
169 
171 {
173  fanCurve_.write(os);
174  os.writeKeyword("direction")
175  << fanFlowDirectionNames_[direction_] << token::END_STATEMENT << nl;
176 }
177 
178 
179 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
180 
181 namespace Foam
182 {
184  (
187  );
188 };
189 
190 
191 // ************************************************************************* //
Foam::surfaceFields.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
surfaceScalarField & phi
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#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
This boundary condition can be applied to assign either a pressure inlet or outlet total pressure con...
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.
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:52
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
makePatchTypeField(fvPatchVectorField, SRFFreestreamVelocityFvPatchVectorField)
const dimensionSet & dimensions() const
Return dimensions.
Type gSum(const FieldField< Field, Type > &f)
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:262
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
const dimensionSet dimDensity
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
static const NamedEnum< fanFlowDirection, 2 > fanFlowDirectionNames_
Fan flow directions names.
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.
fanPressureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
virtual void write(Ostream &) const
Write.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:576
const dimensionSet dimVelocity
This boundary condition provides a total pressure condition. Four variants are possible: ...