All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
geomCellLooper.H
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) 2011-2019 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::geomCellLooper
26 
27 Description
28  Implementation of cellLooper. Does pure geometric cut through cell.
29 
30  Handles all cell shapes in the same way: cut edges with plane through
31  cell centre and normal in direction of provided direction. Snaps cuts
32  close to edge endpoints (close = snapTol * minEdgeLen) to vertices.
33 
34  Currently determines cuts through edges (and edgeendpoints close to plane)
35  in random order and then sorts them acc. to angle. Could be converted to
36  use walk but problem is that face can be cut multiple times (since does
37  not need to be convex). Another problem is that edges parallel to plane
38  might not be cut. So these are handled by looking at the distance from
39  edge endpoints to the plane.
40 
41 SourceFiles
42  geomCellLooper.C
43 
44 \*---------------------------------------------------------------------------*/
45 
46 #ifndef geomCellLooper_H
47 #define geomCellLooper_H
48 
49 #include "cellLooper.H"
50 #include "typeInfo.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 // Forward declaration of classes
58 class plane;
59 
60 /*---------------------------------------------------------------------------*\
61  Class geomCellLooper Declaration
62 \*---------------------------------------------------------------------------*/
63 
64 class geomCellLooper
65 :
66  public cellLooper
67 {
68  // Static
69 
70  //- Tolerance for point equal test. Fraction of edge length.
71  static const scalar pointEqualTol_;
72 
73  //- Tolerance for cut through edges to get snapped to edge end point.
74  // Fraction of length of minimum connected edge length.
75  static scalar snapTol_;
76 
77 
78  // Private Member Functions
79 
80  //- Min length of attached edges
81  scalar minEdgeLen(const label vertI) const;
82 
83  //- Return true and set weight if edge is cut
84  bool cutEdge
85  (
86  const plane& cutPlane,
87  const label edgeI,
88  scalar& weight
89  ) const;
90 
91  //- Snaps cut through edge by cut through vertex (if weight closer than
92  // tol to 0 or 1). Returns vertex label snapped to or -1.
93  label snapToVert
94  (
95  const scalar tol,
96  const label edgeI,
97  const scalar weight
98  ) const;
99 
100  //- Gets two (random) vectors perpendicular to n and each other to be
101  // used as base.
102  void getBase
103  (
104  const vector& n,
105  vector& e0,
106  vector& e1
107  ) const;
108 
109  //- Return true if the cut edge at loop[index] is in between the cuts
110  // through the edge end points.
111  bool edgeEndsCut(const labelList&, const label index) const;
112 
113 
114 public:
115 
116  //- Runtime type information
117  TypeName("geomCellLooper");
118 
119 
120  // Static Functions
122  static scalar snapTol()
123  {
124  return snapTol_;
125  }
127  static void setSnapTol(const scalar tol)
128  {
129  snapTol_ = tol;
130  }
131 
132 
133  // Constructors
134 
135  //- Construct from components
136  geomCellLooper(const polyMesh& mesh);
137 
138  //- Disallow default bitwise copy construction
139  geomCellLooper(const geomCellLooper&) = delete;
140 
141 
142  //- Destructor
143  virtual ~geomCellLooper();
144 
145 
146  // Member Functions
147 
148  //- Create cut along circumference of celli. Gets current mesh cuts.
149  // Cut along circumference is expressed as loop of cuts plus weights
150  // for cuts along edges (only valid for edge cuts).
151  // Return true if successful cut.
152  virtual bool cut
153  (
154  const vector& refDir,
155  const label celli,
156  const boolList& vertIsCut,
157  const boolList& edgeIsCut,
158  const scalarField& edgeWeight,
159 
160  labelList& loop,
161  scalarField& loopWeights
162  ) const;
163 
164  //- Same but now also base point of cut provided (instead of always
165  // cell centre)
166  virtual bool cut
167  (
168  const plane& cutPlane,
169  const label celli,
170  const boolList& vertIsCut,
171  const boolList& edgeIsCut,
172  const scalarField& edgeWeight,
173 
174  labelList& loop,
175  scalarField& loopWeights
176  ) const;
177 
178 
179  // Member Operators
180 
181  //- Disallow default bitwise assignment
182  void operator=(const geomCellLooper&) = delete;
183 };
184 
185 
186 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
187 
188 } // End namespace Foam
189 
190 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
191 
192 #endif
193 
194 // ************************************************************************* //
TypeName("geomCellLooper")
Runtime type information.
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 void setSnapTol(const scalar tol)
Abstract base class. Concrete implementations know how to cut a cell (i.e. determine a loop around th...
Definition: cellLooper.H:69
virtual ~geomCellLooper()
Destructor.
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
Definition: plane.H:60
geomCellLooper(const polyMesh &mesh)
Construct from components.
const polyMesh & mesh() const
Definition: edgeVertex.H:93
void operator=(const geomCellLooper &)=delete
Disallow default bitwise assignment.
virtual bool cut(const vector &refDir, const label celli, const boolList &vertIsCut, const boolList &edgeIsCut, const scalarField &edgeWeight, labelList &loop, scalarField &loopWeights) const
Create cut along circumference of celli. Gets current mesh cuts.
label n
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static scalar snapTol()
Namespace for OpenFOAM.
Implementation of cellLooper. Does pure geometric cut through cell.