smoothDelta.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 | Copyright (C) 2011-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 Class
25  Foam::smoothDelta
26 
27 Description
28  Smoothed delta which takes a given simple geometric delta and applies
29  smoothing to it such that the ratio of deltas between two cells is no
30  larger than a specified amount, typically 1.15.
31 
32 SourceFiles
33  smoothDelta.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef smoothDelta_H
38 #define smoothDelta_H
39 
40 #include "LESdelta.H"
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 namespace LESModels
47 {
48 
49 /*---------------------------------------------------------------------------*\
50  Class smoothDelta Declaration
51 \*---------------------------------------------------------------------------*/
52 
53 class smoothDelta
54 :
55  public LESdelta
56 {
57 public:
58 
59  //- Public member class used by mesh-wave to propagate the delta-ratio
60  class deltaData
61  {
62  scalar delta_;
63 
64  // Private Member Functions
65 
66  //- Update. Gets information from neighbouring face/cell and
67  // uses this to update itself (if necessary) and return true.
68  template<class TrackingData>
69  inline bool update
70  (
71  const deltaData& w2,
72  const scalar scale,
73  const scalar tol,
74  TrackingData& td
75  );
76 
77 
78  public:
79 
80  // Constructors
81 
82  //- Construct null
83  inline deltaData();
84 
85  //- Construct from delta value
86  inline deltaData(const scalar delta);
87 
88 
89  // Member Functions
90 
91  // Access
92 
93  scalar delta() const
94  {
95  return delta_;
96  }
97 
98 
99  // Needed by FaceCellWave
100 
101  //- Check whether origin has been changed at all or
102  // still contains original (invalid) value.
103  template<class TrackingData>
104  inline bool valid(TrackingData& td) const;
105 
106  //- Check for identical geometrical data.
107  // Used for cyclics checking.
108  template<class TrackingData>
109  inline bool sameGeometry
110  (
111  const polyMesh&,
112  const deltaData&,
113  const scalar,
114  TrackingData& td
115  ) const;
116 
117  //- Convert any absolute coordinates into relative to
118  // (patch)face centre
119  template<class TrackingData>
120  inline void leaveDomain
121  (
122  const polyMesh&,
123  const polyPatch&,
124  const label patchFacei,
125  const point& faceCentre,
126  TrackingData& td
127  );
128 
129  //- Reverse of leaveDomain
130  template<class TrackingData>
131  inline void enterDomain
132  (
133  const polyMesh&,
134  const polyPatch&,
135  const label patchFacei,
136  const point& faceCentre,
137  TrackingData& td
138  );
139 
140  //- Apply rotation matrix to any coordinates
141  template<class TrackingData>
142  inline void transform
143  (
144  const polyMesh&,
145  const tensor&,
146  TrackingData& td
147  );
148 
149  //- Influence of neighbouring face.
150  template<class TrackingData>
151  inline bool updateCell
152  (
153  const polyMesh&,
154  const label thisCelli,
155  const label neighbourFacei,
156  const deltaData& neighbourInfo,
157  const scalar tol,
158  TrackingData& td
159  );
160 
161  //- Influence of neighbouring cell.
162  template<class TrackingData>
163  inline bool updateFace
164  (
165  const polyMesh&,
166  const label thisFacei,
167  const label neighbourCelli,
168  const deltaData& neighbourInfo,
169  const scalar tol,
170  TrackingData& td
171  );
172 
173  //- Influence of different value on same face.
174  template<class TrackingData>
175  inline bool updateFace
176  (
177  const polyMesh&,
178  const label thisFacei,
179  const deltaData& neighbourInfo,
180  const scalar tol,
181  TrackingData& td
182  );
183 
184  //- Same (like operator==)
185  template<class TrackingData>
186  inline bool equal(const deltaData&, TrackingData& td) const;
187 
188  // Member Operators
189 
190  // Needed for List IO
191  inline bool operator==(const deltaData&) const;
192 
193  inline bool operator!=(const deltaData&) const;
194 
195  // IOstream Operators
196 
197  friend Ostream& operator<<
198  (
199  Ostream& os,
200  const deltaData& wDist
201  )
202  {
203  return os << wDist.delta_;
204  }
206  friend Istream& operator>>(Istream& is, deltaData& wDist)
207  {
208  return is >> wDist.delta_;
209  }
210  };
211 
212 
213 private:
214 
215  // Private data
216 
217  autoPtr<LESdelta> geometricDelta_;
218  scalar maxDeltaRatio_;
219 
220 
221  // Private Member Functions
222 
223  //- Disallow default bitwise copy construct and assignment
224  smoothDelta(const smoothDelta&);
225  void operator=(const smoothDelta&);
226 
227  // Calculate the delta values
228  void calcDelta();
229 
230  //- Fill changedFaces (with face labels) and changedFacesInfo
231  // (with delta).
232  // This is the initial set of faces from which to start the waves.
233  // Since there might be lots of places with delta jumps we can follow
234  // various strategies for this initial 'seed'.
235  // - start from single cell/face and let FaceCellWave pick up all
236  // others from there. might be quite a few waves before everything
237  // settles.
238  // - start from all faces. Lots of initial transfers.
239  // We do something in between:
240  // - start from all faces where there is a jump. Since we cannot easily
241  // determine this across coupled patches (cyclic, processor)
242  // introduce all faces of these and let FaceCellWave sort it out.
243  void setChangedFaces
244  (
245  const polyMesh& mesh,
246  const volScalarField& delta,
247  DynamicList<label>& changedFaces,
248  DynamicList<deltaData>& changedFacesInfo
249  );
250 
251 
252 public:
253 
254  //- Runtime type information
255  TypeName("smooth");
256 
257 
258  // Constructors
259 
260  //- Construct from name, turbulenceModel and dictionary
262  (
263  const word& name,
265  const dictionary&
266  );
267 
268 
269  //- Destructor
270  virtual ~smoothDelta()
271  {}
272 
273 
274  // Member Functions
275 
276  //- Read the LESdelta dictionary
277  virtual void read(const dictionary&);
278 
279  // Correct values
280  virtual void correct();
281 };
282 
283 
284 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
285 
286 } // End namespace LESModels
287 
288 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
289 
290 //- Data associated with deltaData type are contiguous
291 template<>
292 inline bool contiguous<LESModels::smoothDelta::deltaData>()
293 {
294  return true;
295 }
296 
297 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
298 
299 } // End namespace Foam
300 
301 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
302 
303 #include "smoothDeltaDeltaDataI.H"
304 
305 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
306 
307 #endif
308 
309 // ************************************************************************* //
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
virtual void read(const dictionary &)
Read the LESdelta dictionary.
Definition: smoothDelta.C:170
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
Abstract base class for LES deltas.
Definition: LESdelta.H:50
TypeName("smooth")
Runtime type information.
Abstract base class for turbulence models (RAS, LES and laminar).
void transform(const polyMesh &, const tensor &, TrackingData &td)
Apply rotation matrix to any coordinates.
bool sameGeometry(const polyMesh &, const deltaData &, const scalar, TrackingData &td) const
Check for identical geometrical data.
bool updateFace(const polyMesh &, const label thisFacei, const label neighbourCelli, const deltaData &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring cell.
dynamicFvMesh & mesh
bool operator!=(const deltaData &) const
void leaveDomain(const polyMesh &, const polyPatch &, const label patchFacei, const point &faceCentre, TrackingData &td)
Convert any absolute coordinates into relative to.
bool operator==(const deltaData &) const
A class for handling words, derived from string.
Definition: word.H:59
bool updateCell(const polyMesh &, const label thisCelli, const label neighbourFacei, const deltaData &neighbourInfo, const scalar tol, TrackingData &td)
Influence of neighbouring face.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
friend Istream & operator>>(Istream &is, deltaData &wDist)
Definition: smoothDelta.H:205
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void enterDomain(const polyMesh &, const polyPatch &, const label patchFacei, const point &faceCentre, TrackingData &td)
Reverse of leaveDomain.
bool equal(const deltaData &, TrackingData &td) const
Same (like operator==)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const turbulenceModel & turbulence() const
Return turbulenceModel reference.
Definition: LESdelta.H:131
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual ~smoothDelta()
Destructor.
Definition: smoothDelta.H:269
Namespace for OpenFOAM.
Public member class used by mesh-wave to propagate the delta-ratio.
Definition: smoothDelta.H:59