laplaceFilter.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 "laplaceFilter.H"
29 #include "fvm.H"
30 #include "fvc.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(laplaceFilter, 0);
37  addToRunTimeSelectionTable(LESfilter, laplaceFilter, dictionary);
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 :
45  LESfilter(mesh),
46  widthCoeff_(widthCoeff),
47  coeff_
48  (
49  IOobject
50  (
51  "laplaceFilterCoeff",
52  mesh.time().timeName(),
53  mesh
54  ),
55  mesh,
57  calculatedFvPatchScalarField::typeName
58  )
59 {
60  coeff_.ref() = pow(mesh.V(), 2.0/3.0)/widthCoeff_;
61 }
62 
63 
65 :
66  LESfilter(mesh),
67  widthCoeff_
68  (
69  bd.optionalSubDict(type() + "Coeffs").lookup<scalar>("widthCoeff")
70  ),
71  coeff_
72  (
73  IOobject
74  (
75  "laplaceFilterCoeff",
76  mesh.time().timeName(),
77  mesh
78  ),
79  mesh,
81  calculatedFvPatchScalarField::typeName
82  )
83 {
84  coeff_.ref() = pow(mesh.V(), 2.0/3.0)/widthCoeff_;
85 }
86 
87 
88 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89 
91 {
92  bd.optionalSubDict(type() + "Coeffs").lookup("widthCoeff") >> widthCoeff_;
93 }
94 
95 
96 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
97 
98 Foam::tmp<Foam::volScalarField> Foam::laplaceFilter::operator()
99 (
100  const tmp<volScalarField>& unFilteredField
101 ) const
102 {
103  correctBoundaryConditions(unFilteredField);
104 
105  tmp<volScalarField> filteredField =
106  unFilteredField() + fvc::laplacian(coeff_, unFilteredField());
107 
108  unFilteredField.clear();
109 
110  return filteredField;
111 }
112 
113 
114 Foam::tmp<Foam::volVectorField> Foam::laplaceFilter::operator()
115 (
116  const tmp<volVectorField>& unFilteredField
117 ) const
118 {
119  correctBoundaryConditions(unFilteredField);
120 
121  tmp<volVectorField> filteredField =
122  unFilteredField() + fvc::laplacian(coeff_, unFilteredField());
123 
124  unFilteredField.clear();
125 
126  return filteredField;
127 }
128 
129 
130 Foam::tmp<Foam::volSymmTensorField> Foam::laplaceFilter::operator()
131 (
132  const tmp<volSymmTensorField>& unFilteredField
133 ) const
134 {
135  correctBoundaryConditions(unFilteredField);
136 
137  tmp<volSymmTensorField> filteredField =
138  unFilteredField() + fvc::laplacian(coeff_, unFilteredField());
139 
140  unFilteredField.clear();
141 
142  return filteredField;
143 }
144 
145 
146 Foam::tmp<Foam::volTensorField> Foam::laplaceFilter::operator()
147 (
148  const tmp<volTensorField>& unFilteredField
149 ) const
150 {
151  correctBoundaryConditions(unFilteredField);
152 
153  tmp<volTensorField> filteredField =
154  unFilteredField() + fvc::laplacian(coeff_, unFilteredField());
155 
156  unFilteredField.clear();
157 
158  return filteredField;
159 }
160 
161 
162 // ************************************************************************* //
Abstract class for LES filters.
Definition: LESfilter.H:54
virtual void read(const dictionary &)
Read the LESfilter dictionary.
Definition: laplaceFilter.C:90
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
const fvMesh & mesh() const
Return mesh reference.
Definition: LESfilter.H:130
laplaceFilter(const fvMesh &mesh, scalar widthCoeff)
Construct from components.
Definition: laplaceFilter.C:43
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
Macros for easy insertion into run-time selection tables.
const dimensionSet dimLength
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
const dictionary & optionalSubDict(const word &) const
Find and return a sub-dictionary if found.
Definition: dictionary.C:1040
word timeName
Definition: getTimeIndex.H:3
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
Internal & ref()
Return a reference to the dimensioned internal field.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void correctBoundaryConditions(const tmp< GeoFieldType > &tgf) const
Temporary function to ensure the coupled boundary conditions of the.
Definition: LESfilter.H:72
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
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:92
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844