wallDist.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) 2015-2016 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 "wallDist.H"
27 #include "wallPolyPatch.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33  defineTypeNameAndDebug(wallDist, 0);
34 }
35 
36 
37 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
38 
39 void Foam::wallDist::constructn() const
40 {
41  n_ = tmp<volVectorField>
42  (
43  new volVectorField
44  (
45  IOobject
46  (
47  "n" & patchTypeName_,
48  mesh().time().timeName(),
49  mesh()
50  ),
51  mesh(),
52  dimensionedVector("n" & patchTypeName_, dimless, Zero),
53  patchDistMethod::patchTypes<vector>(mesh(), patchIDs_)
54  )
55  );
56 
57  const fvPatchList& patches = mesh().boundary();
58 
59  volVectorField::Boundary& nbf = n_.ref().boundaryFieldRef();
60 
61  forAllConstIter(labelHashSet, patchIDs_, iter)
62  {
63  label patchi = iter.key();
64  nbf[patchi] == patches[patchi].nf();
65  }
66 }
67 
68 
69 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
70 
71 Foam::wallDist::wallDist(const fvMesh& mesh, const word& patchTypeName)
72 :
74  patchIDs_(mesh.boundaryMesh().findPatchIDs<wallPolyPatch>()),
75  patchTypeName_(patchTypeName),
76  pdm_
77  (
79  (
80  static_cast<const fvSchemes&>(mesh)
81  .subDict(patchTypeName_ & "Dist"),
82  mesh,
83  patchIDs_
84  )
85  ),
86  y_
87  (
88  IOobject
89  (
90  "y" & patchTypeName_,
91  mesh.time().timeName(),
92  mesh
93  ),
94  mesh,
95  dimensionedScalar("y" & patchTypeName_, dimLength, SMALL),
96  patchDistMethod::patchTypes<scalar>(mesh, patchIDs_)
97  ),
98  nRequired_
99  (
100  static_cast<const fvSchemes&>(mesh).subDict(patchTypeName_ & "Dist")
101  .lookupOrDefault<Switch>("nRequired", false)
102  ),
103  n_(volVectorField::null())
104 {
105  if (nRequired_)
106  {
107  constructn();
108  }
109 
110  movePoints();
111 }
112 
113 
114 Foam::wallDist::wallDist
115 (
116  const fvMesh& mesh,
117  const labelHashSet& patchIDs,
118  const word& patchTypeName
119 )
120 :
122  patchIDs_(patchIDs),
123  patchTypeName_(patchTypeName),
124  pdm_
125  (
127  (
128  static_cast<const fvSchemes&>(mesh)
129  .subDict(patchTypeName_ & "Dist"),
130  mesh,
131  patchIDs_
132  )
133  ),
134  y_
135  (
136  IOobject
137  (
138  "y" & patchTypeName_,
139  mesh.time().timeName(),
140  mesh
141  ),
142  mesh,
143  dimensionedScalar("y" & patchTypeName_, dimLength, SMALL),
144  patchDistMethod::patchTypes<scalar>(mesh, patchIDs_)
145  ),
146  nRequired_
147  (
148  static_cast<const fvSchemes&>(mesh).subDict(patchTypeName_ & "Dist")
149  .lookupOrDefault<Switch>("nRequired", false)
150  ),
152 {
153  if (nRequired_)
154  {
155  constructn();
156  }
157 
158  movePoints();
159 }
160 
161 
162 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
163 
165 {}
166 
167 
168 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
169 
171 {
172  if (isNull(n_()))
173  {
175  << "n requested but 'nRequired' not specified in the "
176  << (patchTypeName_ & "Dist") << " dictionary" << nl
177  << " Recalculating y and n fields." << endl;
178 
179  nRequired_ = true;
180  constructn();
181  pdm_->correct(y_, n_.ref());
182  }
183 
184  return n_();
185 }
186 
187 
189 {
190  if (pdm_->movePoints())
191  {
192  if (nRequired_)
193  {
194  return pdm_->correct(y_, n_.ref());
195  }
196  else
197  {
198  return pdm_->correct(y_);
199  }
200  }
201  else
202  {
203  return false;
204  }
205 }
206 
207 
209 {
210  pdm_->updateMesh(mpm);
211  movePoints();
212 }
213 
214 
215 // ************************************************************************* //
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
static autoPtr< patchDistMethod > New(const dictionary &dict, const fvMesh &mesh, const labelHashSet &patchIDs)
wordList patchTypes(nPatches)
virtual ~wallDist()
Destructor.
Definition: wallDist.C:164
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none.
Definition: Switch.H:60
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
HashSet< label, Hash< label > > labelHashSet
A HashSet with label keys.
Definition: HashSet.H:210
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:85
bool isNull(const T &t)
Return true if t is a reference to the nullObject of type T.
Definition: nullObjectI.H:40
dynamicFvMesh & mesh
A class for handling words, derived from string.
Definition: word.H:59
Foam::wallPolyPatch.
Definition: wallPolyPatch.H:48
virtual void updateMesh(const mapPolyMesh &)
Update the y-field when the mesh changes.
Definition: wallDist.C:208
virtual bool movePoints()
Update the y-field when the mesh moves.
Definition: wallDist.C:188
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:91
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:262
defineTypeNameAndDebug(combustionModel, 0)
static const GeometricField< vector, fvPatchField, volMesh > & null()
Return a null geometric field.
const volVectorField & n() const
Return reference to cached normal-to-wall field.
Definition: wallDist.C:170
PtrList< fvPatch > fvPatchList
container classes for fvPatch
Definition: fvPatchList.H:45
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Selector class for finite volume differencing schemes. fvMesh is derived from fvShemes so that all fi...
Definition: fvSchemes.H:50
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Interface to run-time selectable methods to calculate the distance-to-wall and normal-to-wall fields...
Definition: wallDist.H:68
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
Specialisation of patchDist for wall distance calculation.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:545
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243