wallDist.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) 2015-2018 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 {
42  (
43  "n" & patchTypeName_,
44  mesh(),
46  patchDistMethod::patchTypes<vector>(mesh(), patchIDs_)
47  );
48 
49  const fvPatchList& patches = mesh().boundary();
50 
51  volVectorField::Boundary& nbf = n_.ref().boundaryFieldRef();
52 
53  forAllConstIter(labelHashSet, patchIDs_, iter)
54  {
55  label patchi = iter.key();
56  nbf[patchi] == patches[patchi].nf();
57  }
58 }
59 
60 
61 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
62 
63 Foam::wallDist::wallDist(const fvMesh& mesh, const word& patchTypeName)
64 :
66  patchIDs_(mesh.boundaryMesh().findPatchIDs<wallPolyPatch>()),
67  patchTypeName_(patchTypeName),
68  pdm_
69  (
71  (
72  static_cast<const fvSchemes&>(mesh)
73  .subDict(patchTypeName_ & "Dist"),
74  mesh,
75  patchIDs_
76  )
77  ),
78  y_
79  (
80  IOobject
81  (
82  "y" & patchTypeName_,
83  mesh.time().timeName(),
84  mesh
85  ),
86  mesh,
87  dimensionedScalar("y" & patchTypeName_, dimLength, small),
88  patchDistMethod::patchTypes<scalar>(mesh, patchIDs_)
89  ),
90  nRequired_
91  (
92  static_cast<const fvSchemes&>(mesh).subDict(patchTypeName_ & "Dist")
93  .lookupOrDefault<Switch>("nRequired", false)
94  ),
95  n_(volVectorField::null())
96 {
97  if (nRequired_)
98  {
99  constructn();
100  }
101 
102  movePoints();
103 }
104 
105 
107 (
108  const fvMesh& mesh,
109  const labelHashSet& patchIDs,
110  const word& patchTypeName
111 )
112 :
114  patchIDs_(patchIDs),
115  patchTypeName_(patchTypeName),
116  pdm_
117  (
119  (
120  static_cast<const fvSchemes&>(mesh)
121  .subDict(patchTypeName_ & "Dist"),
122  mesh,
123  patchIDs_
124  )
125  ),
126  y_
127  (
128  IOobject
129  (
130  "y" & patchTypeName_,
131  mesh.time().timeName(),
132  mesh
133  ),
134  mesh,
135  dimensionedScalar("y" & patchTypeName_, dimLength, small),
136  patchDistMethod::patchTypes<scalar>(mesh, patchIDs_)
137  ),
138  nRequired_
139  (
140  static_cast<const fvSchemes&>(mesh).subDict(patchTypeName_ & "Dist")
141  .lookupOrDefault<Switch>("nRequired", false)
142  ),
144 {
145  if (nRequired_)
146  {
147  constructn();
148  }
149 
150  movePoints();
151 }
152 
153 
154 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
155 
157 {}
158 
159 
160 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161 
163 {
164  if (isNull(n_()))
165  {
167  << "n requested but 'nRequired' not specified in the "
168  << (patchTypeName_ & "Dist") << " dictionary" << nl
169  << " Recalculating y and n fields." << endl;
170 
171  nRequired_ = true;
172  constructn();
173  pdm_->correct(y_, n_.ref());
174  }
175 
176  return n_();
177 }
178 
179 
181 {
182  if (pdm_->movePoints())
183  {
184  if (nRequired_)
185  {
186  return pdm_->correct(y_, n_.ref());
187  }
188  else
189  {
190  return pdm_->correct(y_);
191  }
192  }
193  else
194  {
195  return false;
196  }
197 }
198 
199 
201 {
202  pdm_->updateMesh(mpm);
203  movePoints();
204 }
205 
206 
207 // ************************************************************************* //
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)
const volVectorField & n() const
Return reference to cached normal-to-wall field.
Definition: wallDist.C:162
wallDist(const fvMesh &mesh, const word &patchTypeName="wall")
Construct from mesh and optional patch type name.
Definition: wallDist.C:63
wordList patchTypes(nPatches)
virtual ~wallDist()
Destructor.
Definition: wallDist.C:156
static tmp< GeometricField< vector, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< vector >> &)
Return a temporary field constructed from name,.
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
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/any.
Definition: Switch.H:60
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
const dimensionSet dimless
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
const dimensionSet dimLength
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:211
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:46
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:200
virtual bool movePoints()
Update the y-field when the mesh moves.
Definition: wallDist.C:180
word timeName
Definition: getTimeIndex.H:3
static const zero Zero
Definition: zero.H:97
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
static const GeometricField< vector, fvPatchField, volMesh > & null()
Return a null geometric field.
PtrList< fvPatch > fvPatchList
container classes for fvPatch
Definition: fvPatchList.H:45
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Selector class for finite volume differencing schemes. fvMesh is derived from fvSchemes so that all f...
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:92
Specialisation of patchDist for wall distance calculation.
Namespace for OpenFOAM.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540