FreeFOAM The Cross-Platform CFD Toolkit
cuttingPlane.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) 1991-2010 OpenCFD Ltd.
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::cuttingPlane
26 
27 Description
28  Constructs plane through mesh.
29 
30  No attempt at resolving degenerate cases. Since the cut faces are
31  usually quite ugly, they will always be triangulated.
32 
33 Note
34  When the cutting plane coincides with a mesh face, the cell edge on the
35  positive side of the plane is taken.
36 
37 SourceFiles
38  cuttingPlane.C
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #ifndef cuttingPlane_H
43 #define cuttingPlane_H
44 
45 #include <OpenFOAM/plane.H>
46 #include <OpenFOAM/pointField.H>
47 #include <OpenFOAM/faceList.H>
48 #include <surfMesh/MeshedSurface.H>
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 
55 class primitiveMesh;
56 
57 /*---------------------------------------------------------------------------*\
58  Class cuttingPlane Declaration
59 \*---------------------------------------------------------------------------*/
60 
62 :
63  public plane,
64  public MeshedSurface<face>
65 {
66  //- Private typedefs for convenience
68 
69  // Private data
70 
71  //- List of cells cut by the plane
72  labelList cutCells_;
73 
74  // Private Member Functions
75 
76  //- Determine cut cells, possibly restricted to a list of cells
77  void calcCutCells
78  (
79  const primitiveMesh&,
80  const scalarField& dotProducts,
81  const UList<label>& cellIdLabels = UList<label>::null()
82  );
83 
84  //- Determine intersection points (cutPoints).
85  void intersectEdges
86  (
87  const primitiveMesh&,
88  const scalarField& dotProducts,
89  List<label>& edgePoint
90  );
91 
92  //- Walk circumference of cell, starting from startEdgeI crossing
93  // only cut edges. Record cutPoint labels in faceVerts.
94  static bool walkCell
95  (
96  const primitiveMesh&,
97  const UList<label>& edgePoint,
98  const label cellI,
99  const label startEdgeI,
100  DynamicList<label>& faceVerts
101  );
102 
103  //- Determine cuts for all cut cells.
104  void walkCellCuts
105  (
106  const primitiveMesh& mesh,
107  const UList<label>& edgePoint
108  );
109 
110 
111 protected:
112 
113  // Constructors
114 
115  //- Construct plane description without cutting
116  cuttingPlane(const plane&);
117 
118  // Protected Member Functions
119 
120  //- recut mesh with existing planeDesc, restricted to a list of cells
121  void reCut
122  (
123  const primitiveMesh&,
124  const UList<label>& cellIdLabels = UList<label>::null()
125  );
126 
127  //- remap action on triangulation or cleanup
128  virtual void remapFaces(const UList<label>& faceMap);
129 
130 public:
131 
132  // Constructors
133 
134  //- Construct from plane and mesh reference,
135  // possibly restricted to a list of cells
137  (
138  const plane&,
139  const primitiveMesh&,
140  const UList<label>& cellIdLabels = UList<label>::null()
141  );
142 
143 
144  // Member Functions
145 
146  //- Return plane used
147  const plane& planeDesc() const
148  {
149  return static_cast<const plane&>(*this);
150  }
151 
152  //- Return List of cells cut by the plane
153  const labelList& cutCells() const
154  {
155  return cutCells_;
156  }
157 
158  //- Return true or false to question: have any cells been cut?
159  bool cut() const
160  {
161  return cutCells_.size();
162  }
163 
164  //- Sample the cell field
165  template<class Type>
166  tmp<Field<Type> > sample(const Field<Type>&) const;
167 
168  template<class Type>
169  tmp<Field<Type> > sample(const tmp<Field<Type> >&) const;
170 
171 
172  // Member Operators
173 
174  void operator=(const cuttingPlane&);
175 };
176 
177 
178 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
179 
180 } // End namespace Foam
181 
182 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
183 
184 #ifdef NoRepository
185 # include "cuttingPlaneTemplates.C"
186 #endif
187 
188 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
189 
190 #endif
191 
192 // ************************ vim: set sw=4 sts=4 et: ************************ //