swirlFanVelocityFvPatchField.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) 2018-2021 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
30 #include "fvPatchFieldMapper.H"
31 #include "volFields.H"
32 #include "surfaceFields.H"
33 #include "unitConversion.H"
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
37 void Foam::swirlFanVelocityFvPatchField::calcFanJump()
38 {
39  if (this->cyclicPatch().owner())
40  {
41  const scalar rpm = rpm_->value(this->db().time().timeOutputValue());
42 
43  const surfaceScalarField& phi =
44  db().lookupObject<surfaceScalarField>(phiName_);
45 
46  const fvPatchField<scalar>& pOwner =
47  patch().lookupPatchField<volScalarField, scalar>(pName_);
48 
49  const label nbrIndex = this->cyclicPatch().neighbPatchID();
50 
51  const fvPatch& nbrPatch = patch().boundaryMesh()[nbrIndex];
52 
53  const fvPatchField<scalar>& pSlave =
54  nbrPatch.lookupPatchField<volScalarField, scalar>(pName_);
55 
56  scalarField deltaP(mag(pOwner - pSlave));
57 
58  if (phi.dimensions() == dimMass/dimTime)
59  {
60  deltaP /=
61  patch().lookupPatchField<volScalarField, scalar>(rhoName_);
62  }
63 
64  const vector axisHat =
65  gSum(patch().nf()*patch().magSf())/gSum(patch().magSf());
66 
67  vectorField tanDir
68  (
69  axisHat ^ (patch().Cf() - origin_)
70  );
71 
72  tanDir /= (mag(tanDir) + SMALL);
73 
74  scalarField magTangU(patch().size(), Zero);
75 
76  if (useRealRadius_)
77  {
78  const vectorField& pCf = patch().Cf();
79 
80  forAll(pCf, i)
81  {
82  const scalar rMag = mag(pCf[i] - origin_);
83 
84  if (rMag > rInner_ && rMag < rOuter_)
85  {
86  magTangU[i] =
87  (
88  deltaP[i]
89  / stabilise
90  (
91  fanEff_ * rMag * rpmToRads(rpm),
92  VSMALL
93  )
94  );
95  }
96  }
97  }
98  else
99  {
100  if (rEff_ <= 0)
101  {
103  << "Effective radius 'rEff' was ill-specified in the "
104  << "dictionary." << nl
105  << exit(FatalError);
106  }
107  magTangU =
108  (
109  deltaP
110  / stabilise
111  (
112  fanEff_ * rEff_ * rpmToRads(rpm),
113  VSMALL
114  )
115  );
116  }
117 
118  // Calculate the tangential velocity
119  const vectorField tangentialVelocity(magTangU*tanDir);
120 
121  this->setJump(tangentialVelocity);
122  }
123 }
124 
125 
126 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
127 
129 (
130  const fvPatch& p,
132 )
133 :
135  phiName_("phi"),
136  pName_("p"),
137  rhoName_("rho"),
138  origin_(),
139  rpm_(nullptr),
140  fanEff_(1),
141  rEff_(0),
142  rInner_(0),
143  rOuter_(0),
144  useRealRadius_(false)
145 {}
146 
147 
149 (
150  const fvPatch& p,
152  const dictionary& dict
153 )
154 :
156  phiName_(dict.getOrDefault<word>("phi", "phi")),
157  pName_(dict.getOrDefault<word>("p", "p")),
158  rhoName_(dict.getOrDefault<word>("rho", "rho")),
159  origin_
160  (
162  (
163  "origin",
164  returnReduce(patch().size(), maxOp<label>())
165  ? gSum(patch().Cf()*patch().magSf())/gSum(patch().magSf())
166  : Zero
167  )
168  ),
169  rpm_
170  (
171  this->cyclicPatch().owner()
172  ? Function1<scalar>::New("rpm", dict, &db())
173  : nullptr
174  ),
175  fanEff_(dict.getOrDefault<scalar>("fanEff", 1)),
176  rEff_(dict.getOrDefault<scalar>("rEff", 0)),
177  rInner_(dict.getOrDefault<scalar>("rInner", 0)),
178  rOuter_(dict.getOrDefault<scalar>("rOuter", 0)),
179  useRealRadius_(dict.getOrDefault("useRealRadius", false))
180 {}
181 
182 
184 (
185  const swirlFanVelocityFvPatchField& rhs,
186  const fvPatch& p,
188  const fvPatchFieldMapper& mapper
189 )
190 :
191  fixedJumpFvPatchField<vector>(rhs, p, iF, mapper),
192  phiName_(rhs.phiName_),
193  pName_(rhs.pName_),
194  rhoName_(rhs.rhoName_),
195  origin_(rhs.origin_),
196  rpm_(rhs.rpm_.clone()),
197  fanEff_(rhs.fanEff_),
198  rEff_(rhs.rEff_),
199  rInner_(rhs.rInner_),
200  rOuter_(rhs.rOuter_),
201  useRealRadius_(rhs.useRealRadius_)
202 {}
203 
204 
206 (
208 )
209 :
211  phiName_(rhs.phiName_),
212  pName_(rhs.pName_),
213  rhoName_(rhs.rhoName_),
214  origin_(rhs.origin_),
215  rpm_(rhs.rpm_.clone()),
216  fanEff_(rhs.fanEff_),
217  rEff_(rhs.rEff_),
218  rInner_(rhs.rInner_),
219  rOuter_(rhs.rOuter_),
220  useRealRadius_(rhs.useRealRadius_)
221 {}
222 
223 
225 (
226  const swirlFanVelocityFvPatchField& rhs,
228 )
229 :
231  phiName_(rhs.phiName_),
232  pName_(rhs.pName_),
233  rhoName_(rhs.rhoName_),
234  origin_(rhs.origin_),
235  rpm_(rhs.rpm_.clone()),
236  fanEff_(rhs.fanEff_),
237  rEff_(rhs.rEff_),
238  rInner_(rhs.rInner_),
239  rOuter_(rhs.rOuter_),
240  useRealRadius_(rhs.useRealRadius_)
241 {}
242 
243 
244 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
245 
247 {
248  if (this->updated())
249  {
250  return;
251  }
252 
253  calcFanJump();
254 }
255 
256 
258 {
260 
261  if (this->cyclicPatch().owner())
262  {
263  os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
264  os.writeEntryIfDifferent<word>("p", "p", pName_);
265  os.writeEntryIfDifferent<word>("rho", "rho", rhoName_);
266  os.writeEntry("origin", origin_);
267 
268  if (rpm_)
269  {
270  rpm_->writeData(os);
271  }
272 
273  os.writeEntryIfDifferent<scalar>("fanEff", 1, fanEff_);
274 
275  if (useRealRadius_)
276  {
277  os.writeEntry("useRealRadius", "true");
278  os.writeEntryIfDifferent<scalar>("rInner", 0, rInner_);
279  os.writeEntryIfDifferent<scalar>("rOuter", 0, rOuter_);
280  }
281  else
282  {
283  os.writeEntryIfDifferent<scalar>("rEff", 0, rEff_);
284  }
285  }
286 }
287 
288 
289 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
290 
291 namespace Foam
292 {
294  (
297  );
298 }
299 
300 // ************************************************************************* //
Foam::fvPatchField< vector >
volFields.H
Foam::maxOp
Definition: ops.H:223
Foam::Ostream::writeEntryIfDifferent
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:248
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::swirlFanVelocityFvPatchField::swirlFanVelocityFvPatchField
swirlFanVelocityFvPatchField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
Definition: swirlFanVelocityFvPatchField.C:129
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::fixedJumpFvPatchField< vector >::setJump
virtual void setJump(const Field< vector > &jump)
Set the jump field.
Definition: fixedJumpFvPatchField.C:147
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::returnReduce
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Definition: PstreamReduceOps.H:94
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
unitConversion.H
Unit conversion functions.
surfaceFields.H
Foam::surfaceFields.
fvPatchFieldMapper.H
Foam::gSum
Type gSum(const FieldField< Field, Type > &f)
Definition: FieldFieldFunctions.C:594
Foam::Function1
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: propellerInfo.H:291
swirlFanVelocityFvPatchField.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::vectorField
Field< vector > vectorField
Specialisation of Field<T> for vector.
Definition: primitiveFieldsFwd.H:54
Foam::fvPatch::boundaryMesh
const fvBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: fvPatch.H:203
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::stabilise
tmp< DimensionedField< scalar, GeoMesh > > stabilise(const DimensionedField< scalar, GeoMesh > &dsf, const dimensioned< scalar > &ds)
Definition: DimensionedScalarField.C:43
Foam::fixedJumpFvPatchField< vector >
Foam::fixedJumpFvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fixedJumpFvPatchField.C:254
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::cyclicFvPatchField< vector >::cyclicPatch
const cyclicFvPatch & cyclicPatch() const
Return local reference cast into the cyclic patch.
Definition: cyclicFvPatchField.H:167
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
Foam::objectRegistry::lookupObject
const Type & lookupObject(const word &name, const bool recursive=false) const
Definition: objectRegistryTemplates.C:434
Foam::swirlFanVelocityFvPatchField::write
virtual void write(Ostream &os) const
Write.
Definition: swirlFanVelocityFvPatchField.C:257
Foam::fvPatch::Cf
const vectorField & Cf() const
Return face centres.
Definition: fvPatch.C:119
Foam::fvPatch::lookupPatchField
const GeometricField::Patch & lookupPatchField(const word &name, const GeometricField *=nullptr, const Type *=nullptr) const
Definition: fvPatchFvMeshTemplates.C:34
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::dimMass
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:51
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::rpmToRads
constexpr scalar rpmToRads(const scalar rpm) noexcept
Conversion from revolutions/minute to radians/sec.
Definition: unitConversion.H:73
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition: surfaceFieldsFwd.H:54
Foam::fvPatchField< vector >::db
const objectRegistry & db() const
Return local objectRegistry.
Definition: fvPatchField.C:206
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::swirlFanVelocityFvPatchField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: swirlFanVelocityFvPatchField.C:246
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::cyclicFvPatch::neighbPatchID
virtual label neighbPatchID() const
Return neighbour.
Definition: cyclicFvPatch.H:100
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::fvPatchField< vector >::patch
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:357
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::swirlFanVelocityFvPatchField
This boundary condition provides a jump condition for U across a cyclic pressure jump condition and a...
Definition: swirlFanVelocityFvPatchField.H:177
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54