atmBoundaryLayer.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) 2014-2016 OpenFOAM Foundation
9  Copyright (C) 2018-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 
29 #include "atmBoundaryLayer.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
39 :
40  initABL_(false),
41  kappa_(0.41),
42  Cmu_(0.09),
43  C1_(0.0),
44  C2_(1.0),
45  ppMin_((boundBox(pp.points())).min()),
46  time_(time),
47  patch_(pp),
48  flowDir_(nullptr),
49  zDir_(nullptr),
50  Uref_(nullptr),
51  Zref_(nullptr),
52  z0_(nullptr),
53  d_(nullptr)
54 {}
55 
56 
58 (
59  const Time& time,
60  const polyPatch& pp,
61  const dictionary& dict
62 )
63 :
64  initABL_(dict.getOrDefault<bool>("initABL", true)),
65  kappa_
66  (
67  dict.getCheckOrDefault<scalar>("kappa", 0.41, scalarMinMax::ge(SMALL))
68  ),
69  Cmu_(dict.getCheckOrDefault<scalar>("Cmu", 0.09, scalarMinMax::ge(SMALL))),
70  C1_(dict.getOrDefault("C1", 0.0)),
71  C2_(dict.getOrDefault("C2", 1.0)),
72  ppMin_((boundBox(pp.points())).min()),
73  time_(time),
74  patch_(pp),
75  flowDir_(Function1<vector>::New("flowDir", dict, &time)),
76  zDir_(Function1<vector>::New("zDir", dict, &time)),
77  Uref_(Function1<scalar>::New("Uref", dict, &time)),
78  Zref_(Function1<scalar>::New("Zref", dict, &time)),
79  z0_(PatchFunction1<scalar>::New(pp, "z0", dict)),
80  d_(PatchFunction1<scalar>::New(pp, "d", dict))
81 {}
82 
83 
85 (
86  const atmBoundaryLayer& abl,
87  const fvPatch& patch,
88  const fvPatchFieldMapper& mapper
89 )
90 :
91  initABL_(abl.initABL_),
92  kappa_(abl.kappa_),
93  Cmu_(abl.Cmu_),
94  C1_(abl.C1_),
95  C2_(abl.C2_),
96  ppMin_(abl.ppMin_),
97  time_(abl.time_),
98  patch_(patch.patch()),
99  flowDir_(abl.flowDir_.clone()),
100  zDir_(abl.zDir_.clone()),
101  Uref_(abl.Uref_.clone()),
102  Zref_(abl.Zref_.clone()),
103  z0_(abl.z0_.clone(patch_)),
104  d_(abl.d_.clone(patch_))
105 {}
106 
107 
109 :
110  initABL_(abl.initABL_),
111  kappa_(abl.kappa_),
112  Cmu_(abl.Cmu_),
113  C1_(abl.C1_),
114  C2_(abl.C2_),
115  ppMin_(abl.ppMin_),
116  time_(abl.time_),
117  patch_(abl.patch_),
118  flowDir_(abl.flowDir_),
119  zDir_(abl.zDir_),
120  Uref_(abl.Uref_),
121  Zref_(abl.Zref_),
122  z0_(abl.z0_.clone(patch_)),
123  d_(abl.d_.clone(patch_))
124 {}
125 
126 
127 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
128 
130 {
131  const scalar t = time_.timeOutputValue();
132  const vector dir(flowDir_->value(t));
133  const scalar magDir = mag(dir);
134 
135  if (magDir < SMALL)
136  {
138  << "magnitude of " << flowDir_->name() << " = " << magDir
139  << " vector must be greater than zero"
140  << abort(FatalError);
141  }
142 
143  return dir/magDir;
144 }
145 
146 
148 {
149  const scalar t = time_.timeOutputValue();
150  const vector dir(zDir_->value(t));
151  const scalar magDir = mag(dir);
152 
153  if (magDir < SMALL)
154  {
156  << "magnitude of " << zDir_->name() << " = " << magDir
157  << " vector must be greater than zero"
158  << abort(FatalError);
159  }
160 
161  return dir/magDir;
162 }
163 
164 
166 {
167  const scalar t = time_.timeOutputValue();
168  const scalar Uref = Uref_->value(t);
169  const scalar Zref = Zref_->value(t);
170 
171  if (Zref < 0)
172  {
174  << "Negative entry in " << Zref_->name() << " = " << Zref
175  << abort(FatalError);
176  }
177 
178  // (derived from RH:Eq. 6, HW:Eq. 5)
179  return kappa_*Uref/log((Zref + z0)/z0);
180 }
181 
182 
184 {
185  if (z0_)
186  {
187  z0_->autoMap(mapper);
188  }
189  if (d_)
190  {
191  d_->autoMap(mapper);
192  }
193 }
194 
195 
197 (
198  const atmBoundaryLayer& abl,
199  const labelList& addr
200 )
201 {
202  if (z0_)
203  {
204  z0_->rmap(abl.z0_(), addr);
205  }
206  if (d_)
207  {
208  d_->rmap(abl.d_(), addr);
209  }
210 }
211 
212 
214 {
215  const scalar t = time_.timeOutputValue();
216  const scalarField d(d_->value(t));
217  const scalarField z0(max(z0_->value(t), ROOTVSMALL));
218  const scalar groundMin = zDir() & ppMin_;
219 
220  // (YGCJ:Table 1, RH:Eq. 6, HW:Eq. 5)
221  scalarField Un
222  (
223  (Ustar(z0)/kappa_)*log(((zDir() & pCf) - groundMin - d + z0)/z0)
224  );
225 
226  return flowDir()*Un;
227 }
228 
229 
231 {
232  const scalar t = time_.timeOutputValue();
233  const scalarField d(d_->value(t));
234  const scalarField z0(max(z0_->value(t), ROOTVSMALL));
235  const scalar groundMin = zDir() & ppMin_;
236 
237  // (YGCJ:Eq. 21; RH:Eq. 7, HW:Eq. 6 when C1=0 and C2=1)
238  return
239  sqr(Ustar(z0))/sqrt(Cmu_)
240  *sqrt(C1_*log(((zDir() & pCf) - groundMin - d + z0)/z0) + C2_);
241 }
242 
243 
245 {
246  const scalar t = time_.timeOutputValue();
247  const scalarField d(d_->value(t));
248  const scalarField z0(max(z0_->value(t), ROOTVSMALL));
249  const scalar groundMin = zDir() & ppMin_;
250 
251  // (YGCJ:Eq. 22; RH:Eq. 8, HW:Eq. 7 when C1=0 and C2=1)
252  return
253  pow3(Ustar(z0))/(kappa_*((zDir() & pCf) - groundMin - d + z0))
254  *sqrt(C1_*log(((zDir() & pCf) - groundMin - d + z0)/z0) + C2_);
255 }
256 
257 
259 {
260  const scalar t = time_.timeOutputValue();
261  const scalarField d(d_->value(t));
262  const scalarField z0(max(z0_->value(t), ROOTVSMALL));
263  const scalar groundMin = zDir() & ppMin_;
264 
265  // (YGJ:Eq. 13)
266  return Ustar(z0)/(kappa_*sqrt(Cmu_)*((zDir() & pCf) - groundMin - d + z0));
267 }
268 
269 
271 {
272  os.writeEntry("initABL", initABL_);
273  os.writeEntry("kappa", kappa_);
274  os.writeEntry("Cmu", Cmu_);
275  os.writeEntry("C1", C1_);
276  os.writeEntry("C2", C2_);
277  flowDir_->writeData(os);
278  zDir_->writeData(os);
279  Uref_->writeData(os);
280  Zref_->writeData(os);
281  if (z0_)
282  {
283  z0_->writeData(os) ;
284  }
285  if (d_)
286  {
287  d_->writeData(os);
288  }
289 }
290 
291 
292 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
293 
294 } // End namespace Foam
295 
296 // ************************************************************************* //
Foam::PrimitivePatch::points
const Field< point_type > & points() const noexcept
Return reference to global points.
Definition: PrimitivePatch.H:299
Foam::TimeState::timeOutputValue
scalar timeOutputValue() const
Return current time value.
Definition: TimeStateI.H:31
Foam::atmBoundaryLayer::Ustar
tmp< scalarField > Ustar(const scalarField &z0) const
Return friction velocity.
Definition: atmBoundaryLayer.C:165
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::atmBoundaryLayer
Base class to set log-law type ground-normal inlet boundary conditions for wind velocity and turbulen...
Definition: atmBoundaryLayer.H:393
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
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
Foam::atmBoundaryLayer::flowDir
vector flowDir() const
Return flow direction.
Definition: atmBoundaryLayer.C:129
Foam::atmBoundaryLayer::k
tmp< scalarField > k(const vectorField &pCf) const
Return the turbulent kinetic energy distribution for the ATM.
Definition: atmBoundaryLayer.C:230
Foam::atmBoundaryLayer::initABL_
bool initABL_
Definition: atmBoundaryLayer.H:401
Foam::atmBoundaryLayer::U
tmp< vectorField > U(const vectorField &pCf) const
Return the velocity distribution for the ATM.
Definition: atmBoundaryLayer.C:213
Foam::atmBoundaryLayer::write
void write(Ostream &) const
Write.
Definition: atmBoundaryLayer.C:270
Foam::Field< scalar >
Foam::pow3
dimensionedScalar pow3(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:89
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::dictionary::getCheckOrDefault
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:209
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::atmBoundaryLayer::omega
tmp< scalarField > omega(const vectorField &pCf) const
Return the specific dissipation rate distribution for the ATM.
Definition: atmBoundaryLayer.C:258
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
Foam::log
dimensionedScalar log(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:262
os
OBJstream os(runTime.globalPath()/outputName)
Foam::PatchFunction1
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
Definition: PatchFunction1.H:60
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::atmBoundaryLayer::atmBoundaryLayer
atmBoundaryLayer(const Time &time, const polyPatch &pp)
Construct null from time.
Definition: atmBoundaryLayer.C:38
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
atmBoundaryLayer.H
Foam::atmBoundaryLayer::epsilon
tmp< scalarField > epsilon(const vectorField &pCf) const
Return the turbulent dissipation rate distribution for the ATM.
Definition: atmBoundaryLayer.C:244
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::MinMax::ge
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
Definition: MinMaxI.H:31
Foam::Vector< scalar >
Foam::atmBoundaryLayer::autoMap
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Definition: atmBoundaryLayer.C:183
Foam::List< label >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::atmBoundaryLayer::zDir
vector zDir() const
Return the ground-normal direction.
Definition: atmBoundaryLayer.C:147
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::atmBoundaryLayer::rmap
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
Definition: atmBoundaryLayer.C:197
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148