anisotropicFilter.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) 2011-2020 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 "anisotropicFilter.H"
29 #include "wallFvPatch.H"
30 #include "fvc.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(anisotropicFilter, 0);
37  addToRunTimeSelectionTable(LESfilter, anisotropicFilter, dictionary);
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const fvMesh& mesh,
46  scalar widthCoeff
47 )
48 :
49  LESfilter(mesh),
50  widthCoeff_(widthCoeff),
51  coeff_
52  (
53  IOobject
54  (
55  "anisotropicFilterCoeff",
56  mesh.time().timeName(),
57  mesh
58  ),
59  mesh,
61  calculatedFvPatchVectorField::typeName
62  )
63 {
64  for (direction d=0; d<vector::nComponents; d++)
65  {
66  coeff_.primitiveFieldRef().replace
67  (
68  d,
69  (1/widthCoeff_)*
70  sqr
71  (
72  2.0*mesh.V()
73  /fvc::surfaceSum(mag(mesh.Sf().component(d)))().primitiveField()
74  )
75  );
76  }
77 }
78 
79 
81 (
82  const fvMesh& mesh,
83  const dictionary& bd
84 )
85 :
86  LESfilter(mesh),
87  widthCoeff_
88  (
89  bd.optionalSubDict(type() + "Coeffs").lookup<scalar>("widthCoeff")
90  ),
91  coeff_
92  (
93  IOobject
94  (
95  "anisotropicFilterCoeff",
96  mesh.time().timeName(),
97  mesh
98  ),
99  mesh,
101  calculatedFvPatchScalarField::typeName
102  )
103 {
104  for (direction d=0; d<vector::nComponents; d++)
105  {
106  coeff_.primitiveFieldRef().replace
107  (
108  d,
109  (1/widthCoeff_)*
110  sqr
111  (
112  2.0*mesh.V()
113  /fvc::surfaceSum(mag(mesh.Sf().component(d)))().primitiveField()
114  )
115  );
116  }
117 }
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
123 {
124  bd.optionalSubDict(type() + "Coeffs").lookup("widthCoeff") >> widthCoeff_;
125 }
126 
127 
128 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
129 
130 Foam::tmp<Foam::volScalarField> Foam::anisotropicFilter::operator()
131 (
132  const tmp<volScalarField>& unFilteredField
133 ) const
134 {
135  correctBoundaryConditions(unFilteredField);
136 
137  tmp<volScalarField> tmpFilteredField =
138  unFilteredField
139  + (
140  coeff_
142  (
143  mesh().Sf()
144  *fvc::snGrad(unFilteredField())
145  )
146  );
147 
148  unFilteredField.clear();
149 
150  return tmpFilteredField;
151 }
152 
153 
154 Foam::tmp<Foam::volVectorField> Foam::anisotropicFilter::operator()
155 (
156  const tmp<volVectorField>& unFilteredField
157 ) const
158 {
159  correctBoundaryConditions(unFilteredField);
160 
161  tmp<volVectorField> tmpFilteredField =
162  unFilteredField
163  + (
164  coeff_
166  (
167  mesh().Sf()
168  *fvc::snGrad(unFilteredField())
169  )
170  );
171 
172  unFilteredField.clear();
173 
174  return tmpFilteredField;
175 }
176 
177 
178 Foam::tmp<Foam::volSymmTensorField> Foam::anisotropicFilter::operator()
179 (
180  const tmp<volSymmTensorField>& unFilteredField
181 ) const
182 {
183  correctBoundaryConditions(unFilteredField);
184 
185  tmp<volSymmTensorField> tmpFilteredField
186  (
188  (
189  "anisotropicFilteredSymmTensorField",
190  mesh(),
191  unFilteredField().dimensions()
192  )
193  );
194 
195  for (direction d=0; d<symmTensor::nComponents; d++)
196  {
197  tmpFilteredField.ref().replace
198  (
199  d, anisotropicFilter::operator()(unFilteredField().component(d))
200  );
201  }
202 
203  unFilteredField.clear();
204 
205  return tmpFilteredField;
206 }
207 
208 
209 Foam::tmp<Foam::volTensorField> Foam::anisotropicFilter::operator()
210 (
211  const tmp<volTensorField>& unFilteredField
212 ) const
213 {
214  correctBoundaryConditions(unFilteredField);
215 
216  tmp<volTensorField> tmpFilteredField
217  (
219  (
220  "anisotropicFilteredTensorField",
221  mesh(),
222  unFilteredField().dimensions()
223  )
224  );
225 
226  for (direction d=0; d<tensor::nComponents; d++)
227  {
228  tmpFilteredField.ref().replace
229  (
230  d, anisotropicFilter::operator()(unFilteredField().component(d))
231  );
232  }
233 
234  unFilteredField.clear();
235 
236  return tmpFilteredField;
237 }
238 
239 
240 // ************************************************************************* //
void replace(const direction, const GeometricField< cmptType, PatchField, GeoMesh > &)
Abstract class for LES filters.
Definition: LESfilter.H:54
const surfaceVectorField & Sf() const
Return cell face area vectors.
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:237
rDeltaT correctBoundaryConditions()
void surfaceIntegrate(Field< Type > &ivf, const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
dimensionedSymmTensor sqr(const dimensionedVector &dv)
uint8_t direction
Definition: direction.H:45
static tmp< GeometricField< symmTensor, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< symmTensor >> &)
Return a temporary field constructed from name,.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
tmp< GeometricField< cmptType, PatchField, GeoMesh > > component(const direction) const
Return a component of the field.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
fvMesh & mesh
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
static const direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:99
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1060
static word timeName(const scalar, const int precision=curPrecision_)
Return time name of given scalar time.
Definition: Time.C:666
static const zero Zero
Definition: zero.H:97
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void read(const dictionary &)
Read the LESfilter dictionary.
tmp< GeometricField< Type, fvPatchField, volMesh > > surfaceSum(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:98
anisotropicFilter(const fvMesh &mesh, scalar widthCoeff)
Construct from components.
void component(FieldField< Field, typename FieldField< Field, Type >::cmptType > &sf, const FieldField< Field, Type > &f, const direction d)
tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > snGrad(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcSnGrad.C:45
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864