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-2015 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_(pTraits<vector>::zero),
38  zDir_(pTraits<vector>::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  (
65  "atmBoundaryLayer(const dictionary&)"
66  ) << "magnitude of n or z must be greater than zero"
67  << abort(FatalError);
68  }
69 
70  // Ensure direction vectors are normalized
71  flowDir_ /= mag(flowDir_);
72  zDir_ /= mag(zDir_);
73 
74  Ustar_ = kappa_*Uref_/(log((Zref_ + z0_)/z0_));
75 }
76 
77 
79 (
80  const atmBoundaryLayer& ptf,
81  const fvPatchFieldMapper& mapper
82 )
83 :
84  flowDir_(ptf.flowDir_),
85  zDir_(ptf.zDir_),
86  kappa_(ptf.kappa_),
87  Cmu_(ptf.Cmu_),
88  Uref_(ptf.Uref_),
89  Zref_(ptf.Zref_),
90  z0_(ptf.z0_, mapper),
91  zGround_(ptf.zGround_, mapper),
92  Ustar_(ptf.Ustar_, mapper)
93 {}
94 
95 
97 :
98  flowDir_(blpvf.flowDir_),
99  zDir_(blpvf.zDir_),
100  kappa_(blpvf.kappa_),
101  Cmu_(blpvf.Cmu_),
102  Uref_(blpvf.Uref_),
103  Zref_(blpvf.Zref_),
104  z0_(blpvf.z0_),
105  zGround_(blpvf.zGround_),
106  Ustar_(blpvf.Ustar_)
107 {}
108 
109 
110 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
111 
113 {
114  z0_.autoMap(m);
115  zGround_.autoMap(m);
116  Ustar_.autoMap(m);
117 }
118 
119 
121 (
122  const atmBoundaryLayer& blptf,
123  const labelList& addr
124 )
125 {
126  z0_.rmap(blptf.z0_, addr);
127  zGround_.rmap(blptf.zGround_, addr);
128  Ustar_.rmap(blptf.Ustar_, addr);
129 }
130 
131 
133 {
134  scalarField Un
135  (
136  (Ustar_/kappa_)
137  *log(((zDir_ & p) - zGround_ + z0_)/z0_)
138  );
139 
140  return flowDir_*Un;
141 }
142 
143 
145 {
146  return sqr(Ustar_)/sqrt(Cmu_);
147 }
148 
149 
151 {
152  return pow3(Ustar_)/(kappa_*((zDir_ & p) - zGround_ + z0_));
153 }
154 
155 
157 {
158  z0_.writeEntry("z0", os) ;
159  os.writeKeyword("flowDir")
160  << flowDir_ << token::END_STATEMENT << nl;
161  os.writeKeyword("zDir")
162  << zDir_ << token::END_STATEMENT << nl;
163  os.writeKeyword("kappa")
164  << kappa_ << token::END_STATEMENT << nl;
165  os.writeKeyword("Cmu")
166  << Cmu_ << token::END_STATEMENT << nl;
167  os.writeKeyword("Uref")
168  << Uref_ << token::END_STATEMENT << nl;
169  os.writeKeyword("Zref")
170  << Zref_ << token::END_STATEMENT << nl;
171  zGround_.writeEntry("zGround", os) ;
172 }
173 
174 
175 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
176 
177 } // End namespace Foam
178 
179 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
tmp< vectorField > U(const vectorField &p) const
Return the velocity distribution for the ATM.
tmp< scalarField > epsilon(const vectorField &p) const
Return the turbulent dissipation rate distribution for the ATM.
dimensionedScalar pow3(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
Foam::fvPatchFieldMapper.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
Namespace for OpenFOAM.
tmp< scalarField > k(const vectorField &p) const
Return the turbulent kinetic energy distribution for the ATM.
dictionary dict
dimensionedScalar log(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:260
stressControl lookup("compactNormalStress") >> compactNormalStress
void autoMap(const FieldMapper &map)
Map from self.
Definition: Field.C:479
volScalarField & p
Definition: createFields.H:51
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:634
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:47
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
void rmap(const UList< Type > &mapF, const labelUList &mapAddressing)
1 to 1 reverse-map from the given field
Definition: Field.C:505
Traits class for primitives.
Definition: pTraits.H:50
error FatalError
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
void rmap(const atmBoundaryLayer &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
atmBoundaryLayer()
Construct null.
void write(Ostream &) const
Write.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
A class for managing temporary objects.
Definition: PtrList.H:118
This class provides functions to evaluate the velocity and turbulence distributions appropriate for a...