general.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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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 "general.H"
31#include "Tuple2.H"
32#include "Switch.H"
33
34// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35
36namespace Foam
37{
38 namespace tabulatedWallFunctions
39 {
42 (
43 tabulatedWallFunction,
44 general,
45 dictionary
46 );
47 }
48}
49
50const Foam::Enum
51<
53>
55({
56 { interpolationType::itLinear, "linear" },
57});
58
59
60// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
61
63{
65
66 // Calculate Reynolds number
67 forAll(uPlus_, i)
68 {
69 Rey[i] = yPlus_[i]*uPlus_[i];
71 {
72 Rey[i] = ::log10(max(ROOTVSMALL, Rey[i]));
73 }
74 }
75
76 // Populate the U+ vs Re table
78 {
79 scalar Re = i*invertedTable_.dx() + invertedTable_.x0();
81 }
82}
83
84
86(
87 const scalar xi,
88 const scalarList& x,
89 const scalarList& fx
90) const
91{
92 switch (interpType_)
93 {
94 case itLinear:
95 {
96 if (xi <= x[0])
97 {
98 return fx[0];
99 }
100 else if (xi >= x.last())
101 {
102 return fx.last();
103 }
104 else
105 {
106 label i2 = 0;
107 while (x[i2] < xi)
108 {
109 i2++;
110 }
111 label i1 = i2 - 1;
112
113 return (xi - x[i1])/(x[i2] - x[i1])*(fx[i2] - fx[i1]) + fx[i1];
114 }
115
116 break;
117 }
118 default:
119 {
121 << "Unknown interpolation method" << nl
122 << abort(FatalError);
123 }
124 }
125
126 return 0.0;
127}
128
129
130// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
131
133(
134 const dictionary& dict,
135 const polyMesh& mesh
136)
137:
138 tabulatedWallFunction(dict, mesh, typeName),
139 interpType_(interpolationTypeNames_.get("interpType", coeffDict_)),
140 log10YPlus_(coeffDict_.get<bool>("log10YPlus")),
141 log10UPlus_(coeffDict_.get<bool>("log10UPlus")),
142 yPlus_(),
143 uPlus_()
144{
145 List<Tuple2<scalar, scalar>> inputTable;
146
147 coeffDict_.readEntry("inputTable", inputTable);
148
149 if (inputTable.size() < 2)
150 {
152 << "Input table must have at least 2 values" << nl
153 << exit(FatalError);
154 }
155
156 yPlus_.setSize(inputTable.size());
157 uPlus_.setSize(inputTable.size());
158
159 forAll(inputTable, i)
160 {
161 if (log10YPlus_)
162 {
163 yPlus_[i] = pow(10, inputTable[i].first());
164 }
165 else
166 {
167 yPlus_[i] = inputTable[i].first();
168 }
169
170 if (log10UPlus_)
171 {
172 uPlus_[i] = pow(10, inputTable[i].second());
173 }
174 else
175 {
176 uPlus_[i] = inputTable[i].second();
177 }
178 }
179
180 invertTable();
181
182 if (debug)
183 {
184 writeData(Info);
185 }
186}
187
188
189// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
190
192{}
193
194
195// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
196
198(
199 const scalar uPlus
200) const
201{
202 return interpolate(uPlus, uPlus_, yPlus_);
203}
204
205
207(
208 const scalar uPlus
209) const
210{
211 return uPlus*yPlus(uPlus);
212}
213
214
215// * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * * //
216
218{
219 if (invertedTable_.log10())
220 {
221 os << "log10(Re), y+, u+:" << endl;
222 forAll(invertedTable_, i)
223 {
224 scalar uPlus = invertedTable_[i];
225 scalar Re = ::log10(this->Re(uPlus));
226 scalar yPlus = this->yPlus(uPlus);
227 os << Re << ", " << yPlus << ", " << uPlus << endl;
228 }
229 }
230 else
231 {
232 os << "Re, y+, u+:" << endl;
233 forAll(invertedTable_, i)
234 {
235 scalar uPlus = invertedTable_[i];
236 scalar Re = this->Re(uPlus);
237 scalar yPlus = this->yPlus(uPlus);
238 os << Re << ", " << yPlus << ", " << uPlus << endl;
239 }
240 }
241}
242
243
244// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
scalar Re() const
Real part of complex number.
Definition: complexI.H:87
Particle-size distribution model wherein random samples are drawn from a given arbitrary probability ...
Definition: general.H:176
Computes the near wall for turbulence models.
Definition: yPlus.H:160
General relative velocity model.
bool interpolate() const noexcept
Same as isPointData()
List< scalar > yPlus_
Input y+ values.
Definition: general.H:112
List< scalar > uPlus_
Input U+ values.
Definition: general.H:115
virtual scalar interpolate(const scalar xi, const scalarList &x, const scalarList &fx) const
Interpolate.
interpolationType
Enumeration listing available interpolation types.
Definition: general.H:91
static const Enum< interpolationType > interpolationTypeNames_
Definition: general.H:95
virtual scalar Re(const scalar uPlus) const
Return Reynolds number as a function of u+.
virtual void invertTable()
Invert the table.
virtual void writeData(Ostream &os) const
Write to Ostream.
uniformInterpolationTable< scalar > invertedTable_
Inverted table.
scalar dx() const
Return the fixed interval.
scalar x0() const
Return the lower limit.
const Switch & log10() const
Return the log10(x) flag.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
bool
Definition: EEqn.H:20
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
scalar uPlus
scalar yPlus
scalar Rey
OBJstream os(runTime.globalPath()/outputName)
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
scalarField Re(const UList< complex > &cf)
Extract real component.
Definition: complexField.C:159
dimensionedScalar log10(const dimensionedScalar &ds)
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:64
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
const bool writeData(pdfDictionary.get< bool >("writeData"))