atmBoundaryLayer.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) 2014-2024 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 
26 #include "atmBoundaryLayer.H"
27 
28 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
29 
30 const Foam::scalar Foam::atmBoundaryLayer::kappaDefault_ = 0.41;
31 
32 const Foam::scalar Foam::atmBoundaryLayer::CmuDefault_ = 0.09;
33 
34 
35 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
36 
37 void Foam::atmBoundaryLayer::init()
38 {
39  if (mag(flowDir_) < small || mag(zDir_) < small)
40  {
42  << "magnitude of n or z must be greater than zero"
43  << abort(FatalError);
44  }
45 
46  // Ensure direction vectors are normalised
47  flowDir_ /= mag(flowDir_);
48  zDir_ /= mag(zDir_);
49 
50  Ustar_ = kappa_*Uref_/(log((Zref_ + z0_)/z0_));
51 }
52 
53 
54 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
55 
57 :
58  flowDir_(Zero),
59  zDir_(Zero),
60  kappa_(0.41),
61  Cmu_(0.09),
62  Uref_(0),
63  Zref_(0),
64  z0_(0),
65  zGround_(0),
66  Ustar_(0),
67  offset_(false),
68  Ulower_(0),
69  kLower_(0),
70  epsilonLower_(0)
71 {}
72 
73 
75 (
76  const vector& flowDir,
77  const vector& zDir,
78  const scalar Uref,
79  const scalar Zref,
80  const scalarField& z0,
81  const scalarField& zGround,
82  const scalar kappa,
83  const scalar Cmu,
84  const scalar Ulower,
85  const scalar kLower,
86  const scalar epsilonLower
87 )
88 :
89  flowDir_(flowDir),
90  zDir_(zDir),
91  kappa_(kappa),
92  Cmu_(Cmu),
93  Uref_(Uref),
94  Zref_(Zref),
95  z0_(z0),
96  zGround_(zGround),
97  Ustar_(z0.size()),
98  offset_(Ulower != 0),
99  Ulower_(Ulower),
100  kLower_(kLower),
101  epsilonLower_(epsilonLower)
102 {
103  init();
104 }
105 
106 
108 (
109  const vectorField& p,
110  const dictionary& dict
111 )
112 :
113  flowDir_(dict.lookup<vector>("flowDir", dimless)),
114  zDir_(dict.lookup<vector>("zDir", dimless)),
115  kappa_(dict.lookupOrDefault<scalar>("kappa", dimless, kappaDefault_)),
116  Cmu_(dict.lookupOrDefault<scalar>("Cmu", dimless, CmuDefault_)),
117  Uref_(dict.lookup<scalar>("Uref", dimVelocity)),
118  Zref_(dict.lookup<scalar>("Zref", dimLength)),
119  z0_("z0", dimLength, dict, p.size()),
120  zGround_("zGround", dimLength, dict, p.size()),
121  Ustar_(p.size()),
122  offset_(dict.found("Ulower")),
123  Ulower_(dict.lookupOrDefault<scalar>("Ulower", dimVelocity, 0)),
124  kLower_(dict.lookupOrDefault<scalar>("kLower", dimEnergy/dimMass, 0)),
125  epsilonLower_
126  (
127  dict.lookupOrDefault<scalar>
128  (
129  "epsilonLower",
131  0
132  )
133  )
134 {
135  init();
136 }
137 
138 
140 (
141  const atmBoundaryLayer& abl,
142  const fieldMapper& mapper
143 )
144 :
145  flowDir_(abl.flowDir_),
146  zDir_(abl.zDir_),
147  kappa_(abl.kappa_),
148  Cmu_(abl.Cmu_),
149  Uref_(abl.Uref_),
150  Zref_(abl.Zref_),
151  z0_(mapper(abl.z0_)),
152  zGround_(mapper(abl.zGround_)),
153  Ustar_(mapper(abl.Ustar_)),
154  offset_(abl.offset_),
155  Ulower_(abl.Ulower_),
156  kLower_(abl.kLower_),
157  epsilonLower_(abl.epsilonLower_)
158 {}
159 
160 
162 :
163  flowDir_(abl.flowDir_),
164  zDir_(abl.zDir_),
165  kappa_(abl.kappa_),
166  Cmu_(abl.Cmu_),
167  Uref_(abl.Uref_),
168  Zref_(abl.Zref_),
169  z0_(abl.z0_),
170  zGround_(abl.zGround_),
171  Ustar_(abl.Ustar_),
172  offset_(abl.offset_),
173  Ulower_(abl.Ulower_),
174  kLower_(abl.kLower_),
175  epsilonLower_(abl.epsilonLower_)
176 {}
177 
178 
179 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
180 
182 (
183  const atmBoundaryLayer& blptf,
184  const fieldMapper& mapper
185 )
186 {
187  mapper(z0_, blptf.z0_);
188  mapper(zGround_, blptf.zGround_);
189  mapper(Ustar_, blptf.Ustar_);
190 }
191 
192 
194 {
195  z0_.reset(blptf.z0_);
196  zGround_.reset(blptf.zGround_);
197  Ustar_.reset(blptf.Ustar_);
198 }
199 
200 
202 (
203  const vectorField& p
204 ) const
205 {
206  const scalarField Un
207  (
208  (Ustar_/kappa_)
209  *log(max((zDir_ & p) - zGround_ + z0_, z0_)/z0_)
210  );
211 
212  if (offset_)
213  {
214  return flowDir_*Un + flowDir_*Ulower_;
215  }
216  else
217  {
218  return flowDir_*Un;
219  }
220 }
221 
222 
224 (
225  const vectorField& p
226 ) const
227 {
229  (
230  sqr(Ustar_)/sqrt(Cmu_)
231  );
232 
233  if (offset_)
234  {
235  const scalarField z((zDir_ & p) - zGround_);
236  tk.ref() = pos0(z)*tk() + neg(z)*kLower_;
237  }
238 
239  return tk;
240 }
241 
242 
244 (
245  const vectorField& p
246 ) const
247 {
248  tmp<scalarField> tepsilon
249  (
250  pow3(Ustar_)/(kappa_*((zDir_ & p) - zGround_ + z0_))
251  );
252 
253  if (offset_)
254  {
255  const scalarField z((zDir_ & p) - zGround_);
256  tepsilon.ref() = pos0(z)*tepsilon() + neg(z)*epsilonLower_;
257  }
258 
259  return tepsilon;
260 }
261 
262 
264 {
265  writeEntry(os, "z0", z0_) ;
266  writeEntry(os, "flowDir", flowDir_);
267  writeEntry(os, "zDir", zDir_);
268  writeEntry(os, "kappa", kappa_);
269  writeEntry(os, "Cmu", Cmu_);
270  writeEntry(os, "Uref", Uref_);
271  writeEntry(os, "Zref", Zref_);
272 
273  if (offset_)
274  {
275  writeEntry(os, "Ulower", Ulower_);
276  writeEntry(os, "kLower", kLower_);
277  writeEntry(os, "epsilonLower", epsilonLower_);
278  }
279 
280  writeEntry(os, "zGround", zGround_);
281 }
282 
283 
284 // ************************************************************************* //
bool found
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
This class provides functions to evaluate the velocity and turbulence distributions appropriate for a...
void map(const atmBoundaryLayer &, const fieldMapper &)
Map the given atmBoundaryLayer onto this atmBoundaryLayer.
void reset(const atmBoundaryLayer &)
Reset the atmBoundaryLayer to the given atmBoundaryLayer.
void write(Ostream &) const
Write.
tmp< vectorField > U(const vectorField &p) const
Return the velocity distribution for the ATM.
atmBoundaryLayer()
Construct null.
tmp< scalarField > k(const vectorField &p) const
Return the turbulent kinetic energy distribution for the ATM.
tmp< scalarField > epsilon(const vectorField &p) const
Return the turbulent dissipation rate distribution for the ATM.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
static const zero Zero
Definition: zero.H:97
dimensionedScalar pos0(const dimensionedScalar &ds)
const dimensionSet dimEnergy
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
Definition: errorManip.H:131
const dimensionSet dimless
const dimensionSet dimLength
dimensionedScalar log(const dimensionedScalar &ds)
const dimensionSet dimTime
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar neg(const dimensionedScalar &ds)
const dimensionSet dimMass
const dimensionSet dimVelocity
error FatalError
dictionary dict
volScalarField & p