motionSmootherAlgo.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::motionSmootherAlgo
26 
27 Description
28  Given a displacement moves the mesh by scaling the displacement back
29  until there are no more mesh errors.
30 
31  Holds displacement field (read upon construction since need boundary
32  conditions) and scaling factor and optional patch number on which to
33  scale back displacement.
34 
35  E.g.
36  \verbatim
37  // Construct iterative mesh mover.
38  motionSmoother meshMover(mesh, labelList(1, patchi));
39 
40  // Set desired displacement:
41  meshMover.displacement() = ..
42 
43  for (label iter = 0; iter < maxIter; iter++)
44  {
45  if (meshMover.scaleMesh(true))
46  {
47  Info<< "Successfully moved mesh" << endl;
48  return true;
49  }
50  }
51  \endverbatim
52 
53 Note
54  - Shared points (parallel): a processor can have points which are part of
55  pp on another processor but have no pp itself (i.e. it has points
56  and/or edges but no faces of pp). Hence we have to be careful when e.g.
57  synchronising displacements that the value from the processor which has
58  faces of pp get priority. This is currently handled in setDisplacement
59  by resetting the internal displacement to zero before doing anything
60  else. The combine operator used will give preference to non-zero
61  values.
62 
63  - Various routines take baffles. These are sets of boundary faces that
64  are treated as a single internal face. This is a hack used to apply
65  movement to internal faces.
66 
67  - Mesh constraints are looked up from the supplied dictionary. (uses
68  recursive lookup)
69 
70 SourceFiles
71  motionSmootherAlgo.C
72  motionSmootherAlgoTemplates.C
73 
74 \*---------------------------------------------------------------------------*/
75 
76 #ifndef motionSmootherAlgo_H
77 #define motionSmootherAlgo_H
78 
79 #include "pointFields.H"
80 #include "HashSet.H"
81 #include "PackedBoolList.H"
82 #include "indirectPrimitivePatch.H"
83 #include "className.H"
84 
85 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
86 
87 namespace Foam
88 {
89 
90 class polyMeshGeometry;
91 class faceSet;
92 
93 /*---------------------------------------------------------------------------*\
94  Class motionSmootherAlgo Declaration
95 \*---------------------------------------------------------------------------*/
96 
98 {
99  // Private class
100 
101  //- To synchronise displacements. We want max displacement since
102  // this is what is specified on pp and internal mesh will have
103  // zero displacement.
104  class maxMagEqOp
105  {
106 
107  public:
108 
109  void operator()(vector& x, const vector& y) const
110  {
111  for (direction i = 0; i < vector::nComponents; i++)
112  {
113  scalar magX = mag(x[i]);
114  scalar magY = mag(y[i]);
115 
116  if (magX < magY)
117  {
118  x[i] = y[i];
119  }
120  else if (magX == magY)
121  {
122  if (y[i] > x[i])
123  {
124  x[i] = y[i];
125  }
126  }
127  }
128  }
129  };
130 
131 
132  // Private data
133 
134  //- Reference to polyMesh. Non-const since we move mesh.
135  polyMesh& mesh_;
136 
137  //- Reference to pointMesh
138  pointMesh& pMesh_;
139 
140  //- Reference to face subset of all adaptPatchIDs
142 
143  //- Displacement field
144  pointVectorField& displacement_;
145 
146  //- Scale factor for displacement
147  pointScalarField& scale_;
148 
149  //- Starting mesh position
150  pointField& oldPoints_;
151 
152 
153  // Internal data
154 
155  //- Indices of fixedValue patches that we're allowed to modify the
156  // displacement on.
157  const labelList adaptPatchIDs_;
158 
159  // Smoothing and checking parameters
160  dictionary paramDict_;
161 
162  //- Is mesh point on boundary or not
163  PackedBoolList isInternalPoint_;
164 
165  //- Is edge master (always except if on coupled boundary and on
166  // lower processor)
167  PackedBoolList isMasterEdge_;
168 
169 
170  // Private Member Functions
171 
172  //- Average of connected points.
173  template<class Type>
175  (
177  const scalarField& edgeWeight
178  ) const;
179 
180  //- Average postion of connected points.
181  template<class Type>
183  (
185  const scalarField& edgeWeight
186  ) const;
187 
188  //- Check constraints
189  template<class Type>
190  static void checkConstraints
191  (
193  );
194 
195  //- Test synchronisation of generic field (not positions!) on points
196  template<class Type, class CombineOp>
197  void testSyncField
198  (
199  const Field<Type>&,
200  const CombineOp& cop,
201  const Type& zero,
202  const scalar maxMag
203  ) const;
204 
205  //- Test synchronisation of points
206  void testSyncPositions(const pointField&, const scalar maxMag) const;
207 
208  static void checkFld(const pointScalarField&);
209 
210  //- Get points used by given faces
211  labelHashSet getPoints(const labelHashSet&) const;
212 
213  //- Calculate per-edge weight
214  tmp<scalarField> calcEdgeWeights(const pointField&) const;
215 
216  //- Explicit smoothing and min on all affected internal points
217  void minSmooth
218  (
219  const scalarField& edgeWeights,
220  const PackedBoolList& isAffectedPoint,
221  const pointScalarField& fld,
222  pointScalarField& newFld
223  ) const;
224 
225  //- Same but only on selected points (usually patch points)
226  void minSmooth
227  (
228  const scalarField& edgeWeights,
229  const PackedBoolList& isAffectedPoint,
230  const labelList& meshPoints,
231  const pointScalarField& fld,
232  pointScalarField& newFld
233  ) const;
234 
235  //- Scale certain (internal) points of a field
236  void scaleField
237  (
238  const labelHashSet& pointLabels,
239  const scalar scale,
241  ) const;
242 
243  //- As above but points have to be in meshPoints as well
244  // (usually to scale patch points)
245  void scaleField
246  (
247  const labelList& meshPoints,
248  const labelHashSet& pointLabels,
249  const scalar scale,
251  ) const;
252 
253  //- Lower on internal points
254  void subtractField
255  (
256  const labelHashSet& pointLabels,
257  const scalar f,
259  ) const;
260 
261  //- As above but points have to be in meshPoints as well
262  // (usually to scale patch points)
263  void subtractField
264  (
265  const labelList& meshPoints,
266  const labelHashSet& pointLabels,
267  const scalar scale,
269  ) const;
270 
271  //- Helper function. Is point internal?
272  bool isInternalPoint(const label pointi) const;
273 
274  //- Given a set of faces that cause smoothing and a number of
275  // iterations determine the maximum set of points who are affected
276  // and the accordingly affected faces.
277  void getAffectedFacesAndPoints
278  (
279  const label nPointIter,
280  const faceSet& wrongFaces,
281 
282  labelList& affectedFaces,
283  PackedBoolList& isAffectedPoint
284  ) const;
285 
286  //- Disallow default bitwise copy construct
288 
289  //- Disallow default bitwise assignment
290  void operator=(const motionSmootherAlgo&);
291 
292 
293 public:
294 
295  ClassName("motionSmootherAlgo");
296 
297  // Constructors
298 
299  //- Construct from mesh, patches to work on and smoothing parameters.
301  (
302  polyMesh&,
303  pointMesh&,
304  indirectPrimitivePatch& pp, // 'outside' points
305  pointVectorField& displacement,
306  pointScalarField& scale,
307  pointField& oldPoints,
308  const labelList& adaptPatchIDs, // patches forming 'outside'
309  const dictionary& paramDict
310  );
311 
312 
313  //- Destructor
315 
316 
317  // Member Functions
318 
319  // Access
320 
321  //- Reference to mesh
322  const polyMesh& mesh() const;
323 
324  //- Reference to pointMesh
325  const pointMesh& pMesh() const;
326 
327  //- Reference to patch
328  const indirectPrimitivePatch& patch() const;
329 
330  //- Patch labels that are being adapted
331  const labelList& adaptPatchIDs() const;
332 
333  const dictionary& paramDict() const;
334 
335 
336 
337  // Edit
338 
339  //- Take over existing mesh position.
340  void correct();
341 
342 
343  //- Set patch fields on patchIDs to be consistent with
344  // all other boundary conditions
345  static void setDisplacementPatchFields
346  (
347  const labelList& patchIDs,
348  pointVectorField& pointDisplacement
349  );
350 
351  //- Set patch fields on displacement to be consistent with
352  // internal values.
354 
355  //- Set displacement field from displacement on patch points.
356  // Modify provided displacement to be consistent with actual
357  // boundary conditions on displacement. Note: resets the
358  // displacement to be 0 on coupled patches beforehand
359  // to make sure shared points
360  // partially on pp (on some processors) and partially not
361  // (on other processors) get the value from pp.
362  static void setDisplacement
363  (
364  const labelList& patchIDs,
365  const indirectPrimitivePatch& pp,
366  pointField& patchDisp,
367  pointVectorField& displacement
368  );
369 
370  void setDisplacement(pointField& patchDisp);
371 
372  //- Special correctBoundaryConditions which evaluates fixedValue
373  // patches first so they get overwritten with any constraint
374  // bc's.
376 
377  //- Apply optional point constraint (2d correction)
378  void modifyMotionPoints(pointField& newPoints) const;
379 
380  //- Get the current points (oldPoints+scale*displacement)
381  tmp<pointField> curPoints() const;
382 
383  //- Set the errorReduction (by how much to scale the displacement
384  // at error locations) parameter. Returns the old value.
385  // Set to 0 (so revert to old mesh) grows out one cell layer
386  // from error faces.
387  scalar setErrorReduction(const scalar);
388 
389  //- Move mesh with given scale. Return true if mesh ok or has
390  // less than nAllow errors, false
391  // otherwise and locally update scale. Smoothmesh=false means only
392  // patch points get moved.
393  // Parallel ok (as long as displacement field is consistent
394  // across patches)
395  bool scaleMesh
396  (
397  labelList& checkFaces,
398  const bool smoothMesh = true,
399  const label nAllow = 0
400  );
401 
402  //- Move mesh (with baffles) with given scale.
403  bool scaleMesh
404  (
405  labelList& checkFaces,
406  const List<labelPair>& baffles,
407  const bool smoothMesh = true,
408  const label nAllow = 0
409  );
410 
411  //- Move mesh with externally provided mesh constraints
412  bool scaleMesh
413  (
414  labelList& checkFaces,
415  const List<labelPair>& baffles,
416  const dictionary& paramDict,
417  const dictionary& meshQualityDict,
418  const bool smoothMesh = true,
419  const label nAllow = 0
420  );
421 
422 
423  //- Update for new mesh geometry
424  void movePoints();
425 
426  //- Update for new mesh topology
427  void updateMesh();
428 
429 
430  //- Check mesh with mesh settings in dict. Collects incorrect faces
431  // in set. Returns true if one or more faces in error.
432  // Parallel ok.
433  static bool checkMesh
434  (
435  const bool report,
436  const polyMesh& mesh,
437  const dictionary& dict,
438  labelHashSet& wrongFaces
439  );
440 
441  //- Check (subset of mesh) with mesh settings in dict.
442  // Collects incorrect faces in set. Returns true if one
443  // or more faces in error. Parallel ok.
444  static bool checkMesh
445  (
446  const bool report,
447  const polyMesh& mesh,
448  const dictionary& dict,
449  const labelList& checkFaces,
450  labelHashSet& wrongFaces
451  );
452 
453  //- Check (subset of mesh including baffles) with mesh settings
454  // in dict. Collects incorrect faces in set. Returns true if one
455  // or more faces in error. Parallel ok.
456  static bool checkMesh
457  (
458  const bool report,
459  const polyMesh& mesh,
460  const dictionary& dict,
461  const labelList& checkFaces,
462  const List<labelPair>& baffles,
463  labelHashSet& wrongFaces
464  );
465 
466  //- Check part of mesh with mesh settings in dict.
467  // Collects incorrect faces in set. Returns true if one or
468  // more faces in error. Parallel ok.
469  static bool checkMesh
470  (
471  const bool report,
472  const dictionary& dict,
473  const polyMeshGeometry&,
474  const labelList& checkFaces,
475  labelHashSet& wrongFaces
476  );
477 
478  //- Check part of mesh including baffles with mesh settings in dict.
479  // Collects incorrect faces in set. Returns true if one or
480  // more faces in error. Parallel ok.
481  static bool checkMesh
482  (
483  const bool report,
484  const dictionary& dict,
485  const polyMeshGeometry&,
486  const labelList& checkFaces,
487  const List<labelPair>& baffles,
488  labelHashSet& wrongFaces
489  );
490 
491  // Helper functions to manipulate displacement vector.
492 
493  //- Fully explicit smoothing of fields (not positions)
494  // of internal points with varying diffusivity.
495  template<class Type>
496  void smooth
497  (
499  const scalarField& edgeWeight,
501  ) const;
502 };
503 
504 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
505 
506 } // End namespace Foam
507 
508 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
509 
510 #ifdef NoRepository
512 #endif
513 
514 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
515 
516 #endif
517 
518 // ************************************************************************* //
static bool checkMesh(const bool report, const polyMesh &mesh, const dictionary &dict, labelHashSet &wrongFaces)
Check mesh with mesh settings in dict. Collects incorrect faces.
dictionary dict
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
uint8_t direction
Definition: direction.H:46
A list of face labels.
Definition: faceSet.H:48
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
static void setDisplacement(const labelList &patchIDs, const indirectPrimitivePatch &pp, pointField &patchDisp, pointVectorField &displacement)
Set displacement field from displacement on patch points.
bool scaleMesh(labelList &checkFaces, const bool smoothMesh=true, const label nAllow=0)
Move mesh with given scale. Return true if mesh ok or has.
labelList pointLabels(nPoints,-1)
Generic GeometricField class.
const indirectPrimitivePatch & patch() const
Reference to patch.
void movePoints()
Update for new mesh geometry.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
void correct()
Take over existing mesh position.
scalar y
A list of faces which address into the list of points.
const labelList & adaptPatchIDs() const
Patch labels that are being adapted.
static const direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:96
scalar setErrorReduction(const scalar)
Set the errorReduction (by how much to scale the displacement.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
void setDisplacementPatchFields()
Set patch fields on displacement to be consistent with.
void smooth(const GeometricField< Type, pointPatchField, pointMesh > &fld, const scalarField &edgeWeight, GeometricField< Type, pointPatchField, pointMesh > &newFld) const
Fully explicit smoothing of fields (not positions)
Updateable mesh geometry and checking routines.
const dictionary & paramDict() const
void modifyMotionPoints(pointField &newPoints) const
Apply optional point constraint (2d correction)
const polyMesh & mesh() const
Reference to mesh.
labelList f(nPoints)
A bit-packed bool list.
ClassName("motionSmootherAlgo")
void correctBoundaryConditions(pointVectorField &) const
Special correctBoundaryConditions which evaluates fixedValue.
Macro definitions for declaring ClassName(), NamespaceName(), etc.
const pointMesh & pMesh() const
Reference to pointMesh.
A class representing the concept of 0 used to avoid unnecessary manipulations for objects that are kn...
Definition: zero.H:49
dimensioned< scalar > mag(const dimensioned< Type > &)
void updateMesh()
Update for new mesh topology.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A class for managing temporary objects.
Definition: PtrList.H:54
Given a displacement moves the mesh by scaling the displacement back until there are no more mesh err...
tmp< pointField > curPoints() const
Get the current points (oldPoints+scale*displacement)
Namespace for OpenFOAM.