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 | Copyright (C) 2014-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 
26 #include "atmBoundaryLayer.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 :
37  flowDir_(Zero),
38  zDir_(Zero),
39  kappa_(0.41),
40  Cmu_(0.09),
41  Uref_(0),
42  Zref_(0),
43  z0_(0),
44  zGround_(0),
45  Ustar_(0)
46 {}
47 
48 
50 :
51  flowDir_(dict.lookup("flowDir")),
52  zDir_(dict.lookup("zDir")),
53  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
54  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
55  Uref_(readScalar(dict.lookup("Uref"))),
56  Zref_(readScalar(dict.lookup("Zref"))),
57  z0_("z0", dict, p.size()),
58  zGround_("zGround", dict, p.size()),
59  Ustar_(p.size())
60 {
61  if (mag(flowDir_) < SMALL || mag(zDir_) < SMALL)
62  {
64  << "magnitude of n or z must be greater than zero"
65  << abort(FatalError);
66  }
67 
68  // Ensure direction vectors are normalized
69  flowDir_ /= mag(flowDir_);
70  zDir_ /= mag(zDir_);
71 
72  Ustar_ = kappa_*Uref_/(log((Zref_ + z0_)/z0_));
73 }
74 
75 
77 (
78  const atmBoundaryLayer& ptf,
79  const fvPatchFieldMapper& mapper
80 )
81 :
82  flowDir_(ptf.flowDir_),
83  zDir_(ptf.zDir_),
84  kappa_(ptf.kappa_),
85  Cmu_(ptf.Cmu_),
86  Uref_(ptf.Uref_),
87  Zref_(ptf.Zref_),
88  z0_(ptf.z0_, mapper),
89  zGround_(ptf.zGround_, mapper),
90  Ustar_(ptf.Ustar_, mapper)
91 {}
92 
93 
95 :
96  flowDir_(blpvf.flowDir_),
97  zDir_(blpvf.zDir_),
98  kappa_(blpvf.kappa_),
99  Cmu_(blpvf.Cmu_),
100  Uref_(blpvf.Uref_),
101  Zref_(blpvf.Zref_),
102  z0_(blpvf.z0_),
103  zGround_(blpvf.zGround_),
104  Ustar_(blpvf.Ustar_)
105 {}
106 
107 
108 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
109 
111 {
112  z0_.autoMap(m);
113  zGround_.autoMap(m);
114  Ustar_.autoMap(m);
115 }
116 
117 
119 (
120  const atmBoundaryLayer& blptf,
121  const labelList& addr
122 )
123 {
124  z0_.rmap(blptf.z0_, addr);
125  zGround_.rmap(blptf.zGround_, addr);
126  Ustar_.rmap(blptf.Ustar_, addr);
127 }
128 
129 
131 {
132  scalarField Un
133  (
134  (Ustar_/kappa_)
135  *log(((zDir_ & p) - zGround_ + z0_)/z0_)
136  );
137 
138  return flowDir_*Un;
139 }
140 
141 
143 {
144  return sqr(Ustar_)/sqrt(Cmu_);
145 }
146 
147 
149 {
150  return pow3(Ustar_)/(kappa_*((zDir_ & p) - zGround_ + z0_));
151 }
152 
153 
155 {
156  z0_.writeEntry("z0", os) ;
157  os.writeKeyword("flowDir")
158  << flowDir_ << token::END_STATEMENT << nl;
159  os.writeKeyword("zDir")
160  << zDir_ << token::END_STATEMENT << nl;
161  os.writeKeyword("kappa")
162  << kappa_ << token::END_STATEMENT << nl;
163  os.writeKeyword("Cmu")
164  << Cmu_ << token::END_STATEMENT << nl;
165  os.writeKeyword("Uref")
166  << Uref_ << token::END_STATEMENT << nl;
167  os.writeKeyword("Zref")
168  << Zref_ << token::END_STATEMENT << nl;
169  zGround_.writeEntry("zGround", os) ;
170 }
171 
172 
173 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
174 
175 } // End namespace Foam
176 
177 // ************************************************************************* //
dictionary dict
This class provides functions to evaluate the velocity and turbulence distributions appropriate for a...
dimensionedScalar log(const dimensionedScalar &ds)
tmp< vectorField > U(const vectorField &p) const
Return the velocity distribution for the ATM.
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar sqrt(const dimensionedScalar &ds)
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:719
void write(Ostream &) const
Write.
stressControl lookup("compactNormalStress") >> compactNormalStress
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:91
tmp< scalarField > epsilon(const vectorField &p) const
Return the turbulent dissipation rate distribution for the ATM.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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
void autoMap(const FieldMapper &map, const bool applyFlip=true)
Map from self.
Definition: Field.C:514
tmp< scalarField > k(const vectorField &p) const
Return the turbulent kinetic energy distribution for the ATM.
dimensionedScalar pow3(const dimensionedScalar &ds)
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:577
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:54
atmBoundaryLayer()
Construct null.
Namespace for OpenFOAM.