SRFFreestreamVelocityFvPatchVectorField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2021 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 #include "SRFModel.H"
31 #include "steadyStateDdtScheme.H"
32 
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
38 (
39  const fvPatch& p,
41 )
42 :
43  inletOutletFvPatchVectorField(p, iF),
44  relative_(false),
45  UInf_(Zero)
46 {}
47 
48 
51 (
52  const fvPatch& p,
54  const dictionary& dict
55 )
56 :
57  inletOutletFvPatchVectorField(p, iF),
58  relative_(dict.lookupOrDefault("relative", false)),
59  UInf_(dict.lookup("UInf"))
60 {
61  this->phiName_ = dict.lookupOrDefault<word>("phi","phi");
62 
63  fvPatchVectorField::operator=(vectorField("value", dict, p.size()));
64 }
65 
66 
69 (
71  const fvPatch& p,
73  const fvPatchFieldMapper& mapper
74 )
75 :
76  inletOutletFvPatchVectorField(ptf, p, iF, mapper),
77  relative_(ptf.relative_),
78  UInf_(ptf.UInf_)
79 {}
80 
81 
84 (
87 )
88 :
89  inletOutletFvPatchVectorField(srfvpvf, iF),
90  relative_(srfvpvf.relative_),
91  UInf_(srfvpvf.UInf_)
92 {}
93 
94 
95 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
96 
98 {
99  if (updated())
100  {
101  return;
102  }
103 
104  // Get reference to the SRF model
105  const SRF::SRFModel& srf =
106  db().lookupObject<SRF::SRFModel>("SRFProperties");
107 
108 
109  word ddtScheme
110  (
111  this->internalField().mesh()
112  .ddtScheme(this->internalField().name())
113  );
114 
116  {
117  // If not relative to the SRF include the effect of the SRF
118  if (!relative_)
119  {
120  refValue() = UInf_ - srf.velocity(patch().Cf());
121  }
122  // If already relative to the SRF simply supply the inlet value
123  // as a fixed value
124  else
125  {
126  refValue() = UInf_;
127  }
128  }
129  else
130  {
131  scalar time = this->db().time().value();
132  scalar theta = time*mag(srf.omega().value());
133 
134  refValue() =
135  cos(theta)*UInf_ + sin(theta)*(srf.axis() ^ UInf_)
136  - srf.velocity(patch().Cf());
137  }
138 
139  // Set the inlet-outlet choice based on the direction of the freestream
140  valueFraction() = 1.0 - pos0(refValue() & patch().Sf());
141 
143 }
144 
145 
147 {
149  writeEntry(os, "relative", relative_);
150  writeEntry(os, "UInf", UInf_);
151  writeEntry(os, "phi", this->phiName_);
152  writeEntry(os, "value", *this);
153 }
154 
155 
156 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
157 
158 namespace Foam
159 {
161  (
164  );
165 }
166 
167 
168 // ************************************************************************* //
Foam::surfaceFields.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
volVectorField vectorField(fieldObject, mesh)
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:260
Macros for easy insertion into run-time selection tables.
const vector & axis() const
Return the axis of rotation.
Definition: SRFModel.C:109
dynamicFvMesh & mesh
dimensionedScalar cos(const dimensionedScalar &ds)
A class for handling words, derived from string.
Definition: word.H:59
Foam::fvPatchFieldMapper.
SRFFreestreamVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
const Type & value() const
Return const reference to value.
const dimensionedVector & omega() const
Return the angular velocity field [rad/s].
Definition: SRFModel.C:115
static const zero Zero
Definition: zero.H:97
SteadyState implicit/explicit ddt which returns 0.
virtual label size() const
Return size.
Definition: fvPatch.H:156
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedScalar sin(const dimensionedScalar &ds)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:209
Freestream velocity condition to be used in conjunction with the single rotating frame (SRF) model (s...
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Top level model for single rotating frame.
Definition: SRFModel.H:62
vectorField velocity(const vectorField &positions) const
Return velocity vector from positions.
Definition: SRFModel.C:151
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844