twoDPointCorrector.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-2015 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::twoDPointCorrector
26 
27 Description
28  Class applies a two-dimensional correction to mesh motion point field.
29 
30  The correction guarantees that the mesh does not get twisted during motion
31  and thus introduce a third dimension into a 2-D problem.
32 
33  The operation is performed by looping through all edges approximately
34  normal to the plane and enforcing their orthogonality onto the plane by
35  adjusting points on their ends.
36 
37 SourceFiles
38  twoDPointCorrector.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef twoDPointCorrector_H
43 #define twoDPointCorrector_H
44 
45 #include "MeshObject.H"
46 #include "pointField.H"
47 #include "labelList.H"
48 #include "vector.H"
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 
55 // Forward class declarations
56 class polyMesh;
57 
58 /*---------------------------------------------------------------------------*\
59  Class twoDPointCorrector Declaration
60 \*---------------------------------------------------------------------------*/
61 
63 :
64  public MeshObject<polyMesh, UpdateableMeshObject, twoDPointCorrector>
65 {
66  // Private data
67 
68  //- Is 2D correction required, i.e. is the mesh
69  bool required_;
70 
71  //- 2-D plane unit normal
72  mutable vector* planeNormalPtr_;
73 
74  //- Indices of edges normal to plane
75  mutable labelList* normalEdgeIndicesPtr_;
76 
77  //- Flag to indicate a wedge geometry
78  mutable bool isWedge_;
79 
80  //- Wedge axis (if wedge geometry)
81  mutable vector wedgeAxis_;
82 
83  //- Wedge angle (if wedge geometry)
84  mutable scalar wedgeAngle_;
85 
86 
87  // Private Member Functions
88 
89  //- Disallow default bitwise copy construct
91 
92  //- Disallow default bitwise assignment
93  void operator=(const twoDPointCorrector&);
94 
95 
96  //- Calculate addressing
97  void calcAddressing() const;
98 
99  //- Clear addressing
100  void clearAddressing() const;
101 
102  //- Snap a point to the wedge patch(es)
103  void snapToWedge(const vector& n, const point& A, point& p) const;
104 
105 
106  // Static data members
107 
108  //- Edge orthogonality tolerance
109  static const scalar edgeOrthogonalityTol;
110 
111 
112 public:
113 
114  // Declare name of the class and its debug switch
115  ClassName("twoDPointCorrector");
116 
117 
118  // Constructors
119 
120  //- Construct from components
122 
123 
124  //- Destructor
126 
127 
128  // Member Functions
129 
130  //- Is 2D correction required, i.e. is the mesh a wedge or slab
131  bool required() const
132  {
133  return required_;
134  }
135 
136  //- Return plane normal
137  const vector& planeNormal() const;
138 
139  //- Return indices of normal edges.
140  const labelList& normalEdgeIndices() const;
141 
142  //- Return direction normal to plane
143  direction normalDir() const;
144 
145  //- Correct motion points
146  void correctPoints(pointField& p) const;
147 
148  //- Correct motion displacements
149  void correctDisplacement(const pointField& p, vectorField& disp) const;
150 
151  //- Update topology
152  void updateMesh(const mapPolyMesh&);
153 
154  //- Correct weighting factors for moving mesh.
155  bool movePoints();
156 };
157 
158 
159 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
160 
161 } // End namespace Foam
162 
163 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
164 
165 #endif
166 
167 // ************************************************************************* //
void updateMesh(const mapPolyMesh &)
Update topology.
void correctDisplacement(const pointField &p, vectorField &disp) const
Correct motion displacements.
uint8_t direction
Definition: direction.H:46
const vector & planeNormal() const
Return plane normal.
ClassName("twoDPointCorrector")
bool required() const
Is 2D correction required, i.e. is the mesh a wedge or slab.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:85
Class applies a two-dimensional correction to mesh motion point field.
const labelList & normalEdgeIndices() const
Return indices of normal edges.
bool movePoints()
Correct weighting factors for moving mesh.
label n
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
direction normalDir() const
Return direction normal to plane.
volScalarField & p
void correctPoints(pointField &p) const
Correct motion points.
Namespace for OpenFOAM.