turbulentDigitalFilterInletFvPatchField.H
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) 2019-2022 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::turbulentDigitalFilterInletFvPatchField
28
29Group
30 grpInletBoundaryConditions
31
32Description
33 Digital-filter based boundary condition for vector- and scalar-based
34 quantities (e.g. \c U or \c T) to generate synthetic turbulence-alike
35 time-series from input turbulence statistics for LES and DES
36 turbulent flow computations.
37
38 References:
39 \verbatim
40 Digital-filter method-based generator (DFM) (tag:KSJ):
41 Klein, M., Sadiki, A., & Janicka, J. (2003).
42 A digital filter based generation of inflow data for spatially
43 developing direct numerical or large eddy simulations.
44 Journal of computational Physics, 186(2), 652-665.
45 DOI:10.1016/S0021-9991(03)00090-1
46
47 Forward-stepwise method-based generator (FSM) (tag:XC)
48 Xie, Z. T., & Castro, I. P. (2008).
49 Efficient generation of inflow conditions for
50 large eddy simulation of street-scale flows.
51 Flow, turbulence and combustion, 81(3), 449-470.
52 DOI:10.1007/s10494-008-9151-5
53
54 Mass-inflow rate correction (tag:KCX):
55 Kim, Y., Castro, I. P., & Xie, Z. T. (2013).
56 Divergence-free turbulence inflow conditions for
57 large-eddy simulations with incompressible flow solvers.
58 Computers & Fluids, 84, 56-68.
59 DOI:10.1016/j.compfluid.2013.06.001
60
61 Scalar fluctuations (tag:XHW):
62 Xie, Z. T., Hayden, P., & Wood, C. R. (2013).
63 Large-eddy simulation of approaching-flow stratification
64 on dispersion over arrays of buildings.
65 Atmospheric Environment, 71, 64-74.
66 DOI:10.1016/j.atmosenv.2013.01.054
67
68 Okaze, T., & Mochida, A. (2017).
69 Cholesky decomposition–based generation of artificial
70 inflow turbulence including scalar fluctuation.
71 Computers & Fluids, 159, 23-32.
72 DOI:10.1016/j.compfluid.2017.09.005
73 \endverbatim
74
75 In \c DFM or \c FSM, a random number set (mostly white noise), and a group
76 of target statistics (mostly mean flow, Reynolds-stress tensor profiles and
77 integral-scale sets) are merged into a new number set (stochastic
78 time-series, yet consisting of the statistics) by a chain of mathematical
79 operations whose characteristics are designated by the target statistics,
80 so that the realised statistics of the new sets could match the target.
81
82 \verbatim
83 Random number sets ---->-|
84 |
85 DFM or FSM ---> New stochastic time-series consisting
86 | turbulence statistics
87 Turbulence statistics ->-|
88 \endverbatim
89
90 The main difference between \c DFM and \c FSM is that \c FSM replaces
91 the expensive-to-run streamwise convolution summation in \c DFM by a simpler
92 and an almost-equivalent-in-effect numerical procedure in order to reduce
93 computational costs. Accordingly, \c FSM potentially brings computational
94 resource advantages for computations involving relatively large streamwise
95 integral-scale sets and small time-steps.
96
97 Synthetic turbulence is generated on a virtual rectangular structured-mesh
98 plane and is mapped onto this patch by the selected mapping method.
99
100Usage
101 Example of the boundary condition specification:
102 \verbatim
103 <patchName>
104 {
105 // Mandatory entries
106 type turbulentDigitalFilterInlet;
107 n (<label> <label>);
108 mean <PatchFunction1<vector>> or <PatchFunction1<scalar>>;
109 R <PatchFunction1<symmTensor>> or <PatchFunction1<scalar>>;
110 L <tensor> or <vector>;
111
112 // Optional entries
113 kernelType <word>;
114 AMIMethod <word>;
115 coordinateSystem <coordinateSystem>;
116 fsm <bool>;
117
118 // Inherited entries
119 ...
120 }
121 \endverbatim
122
123 where the entries mean:
124 \table
125 Property | Description | Type | Reqd | Deflt
126 type | Type name: turbulentDigitalFilterInlet | word | yes | -
127 n | Number of cells on turbulence generation plane <!--
128 --> (width height) or (e2 e3) | Vector2D<label> | yes | -
129 mean | Mean quantity | PatchFunction1<vector> <!--
130 --> or PatchFunction1<scalar> | yes | -
131 R | Reynolds stresses | PatchFunction1<symmTensor> <!--
132 --> or PatchFunction1<scalar> | yes | -
133 L | Integral scales | tensor or vector | yes | -
134 kernelType | Autocorrelation function form | word | no | Gaussian
135 AMIMethod | Mapping method | word | no | faceAreaWeightAMI
136 coordinateSystem | Coordinate system of turbulence generation plane <!--
137 --> | coordinateSystem | no | -
138 fsm | Flag to turn on the forward-stepwise method | bool | no | false
139 \endtable
140
141 The inherited entries are elaborated in:
142 - \link fixedValueFvPatchFields.H \endlink
143 - \link PatchFunction1.H \endlink
144 - \link AMIPatchToPatchInterpolation.H \endlink
145 - \link coordinateSystem.H \endlink
146 - \link IntegralScaleBox.H \endlink
147
148 Options for the \c kernelType entry:
149 \verbatim
150 Gaussian | Standard Gaussian function
151 exponential | Exponential function
152 \endverbatim
153
154Note
155 - The order of elements of the Reynolds stresses for the vector-based
156 condition is "(Rxx Rxy Rxz Ryy Ryz Rzz)" and for the scalar-based
157 condition is "R".
158 - The order of elements of the integral scales for the vector-based
159 condition is "(Lxu Lxv Lxw Lyu Lyv Lyw Lzu Lzv Lzw)" and for the
160 scalar-based condition is "(Lxi Lyi Lzi)", where \c x here indicates
161 the streamwise components and \c y/z lateral components. The streamwise
162 components (i.e. \c x) are Lagrangian time scales in units of seconds
163 and the remaining components are length scales in units of meters.
164 - \c adjustTimeStep=true option is only supported by
165 the forward-stepwise method.
166 - The default global Cartesian coordinate system "(e1 e2 e3)"
167 corresponds to "(x y z)".
168
169See also
170 - turbulentDFSEMInletFvPatchVectorField.C
171
172SourceFiles
173 turbulentDigitalFilterInletFvPatchField.C
174
175\*---------------------------------------------------------------------------*/
176
177#ifndef Foam_turbulentDigitalFilterInletFvPatchField_H
178#define Foam_turbulentDigitalFilterInletFvPatchField_H
179
181#include "PatchFunction1.H"
182#include "IntegralScaleBox.H"
184
185// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
186
187namespace Foam
188{
189
190/*---------------------------------------------------------------------------*\
191 Class turbulentDigitalFilterInletFvPatchField Declaration
192\*---------------------------------------------------------------------------*/
193
194template<class Type>
195class turbulentDigitalFilterInletFvPatchField
196:
197 public fixedValueFvPatchField<Type>
198{
199 // Private Typedefs
200
201 //- Type of input set of Reynolds stresses
202 typedef typename std::conditional
203 <
204 std::is_same<Type, scalar>::value,
205 scalar,
207 >::type TypeR;
208
209
210 // Private Data
211
212 //- Pointer to AMI interpolator
213 mutable autoPtr<AMIPatchToPatchInterpolation> AMIPtr_;
214
215 //- Mean field on patch in global coordinates
216 autoPtr<PatchFunction1<Type>> meanPtr_;
217
218 //- Reynolds stresses field on patch in global coordinates
219 autoPtr<PatchFunction1<TypeR>> Rptr_;
220
221 //- Current time index
222 label curTimeIndex_;
223
224 //- Patch normal into the domain
225 vector patchNormal_;
226
227 //- Integral-scale field container in global coordinates
228 turbulence::IntegralScaleBox<Type> L_;
229
230
231 // Private Member Functions
232
233 //- Return the patch normal
234 vector calcPatchNormal() const;
235
236 //- Initialise the patch content
237 void initialisePatch();
238
239 //- Map two-point correlations (integral scales)
240 void mapL(Field<Type>& fld);
241
242 //- Map one-point correlations (Reynolds stresses) for scalar fields
243 void mapR(scalarField& fld) const;
244
245 //- Map one-point correlations (Reynolds stresses) for vector fields
246 void mapR(vectorField& fld) const;
247
248 //- Map mean field for scalar fields
249 void mapMean(scalarField& fld) const;
250
251 //- Map mean field for vector fields
252 void mapMean(vectorField& fld) const;
253
254
255public:
256
257 //- Runtime type information
258 TypeName("turbulentDigitalFilterInlet");
259
260
261 // Constructors
263 //- Construct from patch and internal field
265 (
266 const fvPatch&,
268 );
269
270 //- Construct from patch, internal field and patch field
272 (
273 const fvPatch&,
275 const Field<Type>&
276 );
277
278 //- Construct from patch, internal field and dictionary
280 (
281 const fvPatch&,
283 const dictionary&
284 );
285
286 //- Construct by mapping given
287 //- turbulentDigitalFilterInletFvPatchField onto a new patch
289 (
291 const fvPatch&,
293 const fvPatchFieldMapper&
294 );
295
296 //- Construct as copy
298 (
300 );
301
302 //- Construct and return a clone
303 virtual tmp<fvPatchField<Type>> clone() const
304 {
306 (
308 );
309 }
310
311 //- Construct as copy setting internal field reference
313 (
316 );
317
318 //- Construct and return a clone setting internal field reference
320 (
322 ) const
323 {
327 );
328 }
329
330
331 // Member Functions
332
333 // Mapping
334
335 //- Map (and resize as needed) from self given a mapping object
336 virtual void autoMap(const fvPatchFieldMapper& m);
337
338 //- Reverse map the given fvPatchField onto this fvPatchField
339 virtual void rmap
340 (
341 const fvPatchField<Type>& ptf,
342 const labelList& addr
343 );
344
345
346 // Evaluation
347
348 //- Update the coefficients associated with the patch field
349 virtual void updateCoeffs();
350
351
352 // I-O
353
354 //- Write boundary condition settings
355 virtual void write(Ostream&) const;
356};
357
358// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
359
360} // End namespace Foam
361
362// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
363
364#ifdef NoRepository
366#endif
367
368// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
369
370#endif
371
372// ************************************************************************* //
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Generic templated field type.
Definition: Field.H:82
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
This boundary condition supplies a fixed value constraint, and is the base class for a number of othe...
A FieldMapper for finite-volume patch fields.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:82
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:71
A class for managing temporary objects.
Definition: tmp.H:65
Class to describe the integral-scale container being used in the turbulentDigitalFilterInletFvPatchFi...
Digital-filter based boundary condition for vector- and scalar-based quantities (e....
TypeName("turbulentDigitalFilterInlet")
Runtime type information.
turbulentDigitalFilterInletFvPatchField(const fvPatch &, const DimensionedField< Type, volMesh > &, const Field< Type > &)
Construct from patch, internal field and patch field.
virtual void rmap(const fvPatchField< Type > &ptf, const labelList &addr)
Reverse map the given fvPatchField onto this fvPatchField.
virtual void autoMap(const fvPatchFieldMapper &m)
Map (and resize as needed) from self given a mapping object.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual tmp< fvPatchField< Type > > clone() const
Construct and return a clone.
virtual tmp< fvPatchField< Type > > clone(const DimensionedField< Type, volMesh > &iF) const
Construct and return a clone setting internal field reference.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
Namespace for OpenFOAM.
SymmTensor< scalar > symmTensor
SymmTensor of scalars, i.e. SymmTensor<scalar>.
Definition: symmTensor.H:59
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Field< vector > vectorField
Specialisation of Field<T> for vector.
runTime write()
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73