FreeFOAM The Cross-Platform CFD Toolkit
isoSurfaceCell.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::isoSurfaceCell
26 
27 Description
28  A surface formed by the iso value.
29  After "Polygonising A Scalar Field Using Tetrahedrons", Paul Bourke
30  and
31  "Regularised Marching Tetrahedra: improved iso-surface extraction",
32  G.M. Treece, R.W. Prager and A.H. Gee.
33 
34  See isoSurface. This is a variant which does tetrahedrisation from
35  triangulation of face to cell centre instead of edge of face to two
36  neighbouring cell centres. This gives much lower quality triangles
37  but they are local to a cell.
38 
39 SourceFiles
40  isoSurfaceCell.C
41 
42 \*---------------------------------------------------------------------------*/
43 
44 #ifndef isoSurfaceCell_H
45 #define isoSurfaceCell_H
46 
47 #include <triSurface/triSurface.H>
48 #include <OpenFOAM/labelPair.H>
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 
57 class polyMesh;
58 
59 /*---------------------------------------------------------------------------*\
60  Class isoSurfaceCell Declaration
61 \*---------------------------------------------------------------------------*/
62 
64 :
65  public triSurface
66 {
67  // Private data
68 
69  enum segmentCutType
70  {
71  NEARFIRST, // intersection close to e.first()
72  NEARSECOND, // ,, e.second()
73  NOTNEAR // no intersection
74  };
75 
76  enum cellCutType
77  {
78  NOTCUT, // not cut
79  SPHERE, // all edges to cell centre cut
80  CUT // normal cut
81  };
82 
83 
84  //- Reference to mesh
85  const polyMesh& mesh_;
86 
87  //- isoSurfaceCell value
88  const scalar iso_;
89 
90  //- When to merge points
91  const scalar mergeDistance_;
92 
93  //- Whether cell might be cut
94  List<cellCutType> cellCutType_;
95 
96  //- Estimated number of cut cells
97  label nCutCells_;
98 
99  //- For every triangle the original cell in mesh
100  labelList meshCells_;
101 
102  //- For every unmerged triangle point the point in the triSurface
103  labelList triPointMergeMap_;
104 
105 
106  // Private Member Functions
107 
108  //- Get location of iso value as fraction inbetween s0,s1
109  scalar isoFraction
110  (
111  const scalar s0,
112  const scalar s1
113  ) const;
114 
115  //List<triFace> triangulate(const face& f) const;
116 
117  //- Determine whether cell is cut
118  cellCutType calcCutType
119  (
120  const PackedBoolList&,
121  const scalarField& cellValues,
122  const scalarField& pointValues,
123  const label
124  ) const;
125 
126  void calcCutTypes
127  (
128  const PackedBoolList&,
129  const scalarField& cellValues,
130  const scalarField& pointValues
131  );
132 
133  static labelPair findCommonPoints
134  (
135  const labelledTri&,
136  const labelledTri&
137  );
138 
139  static point calcCentre(const triSurface&);
140 
141  static pointIndexHit collapseSurface
142  (
145  );
146 
147  //- Determine per cc whether all near cuts can be snapped to single
148  // point.
149  void calcSnappedCc
150  (
151  const PackedBoolList& isTet,
152  const scalarField& cVals,
153  const scalarField& pVals,
154  DynamicList<point>& triPoints,
155  labelList& snappedCc
156  ) const;
157 
158  //- Generate triangles for face connected to pointI
159  void genPointTris
160  (
161  const scalarField& cellValues,
162  const scalarField& pointValues,
163  const label pointI,
164  const face& f,
165  const label cellI,
166  DynamicList<point, 64>& localTriPoints
167  ) const;
168 
169  //- Generate triangles for tet connected to pointI
170  void genPointTris
171  (
172  const scalarField& pointValues,
173  const label pointI,
174  const label cellI,
175  DynamicList<point, 64>& localTriPoints
176  ) const;
177 
178  //- Determine per point whether all near cuts can be snapped to single
179  // point.
180  void calcSnappedPoint
181  (
182  const PackedBoolList& isBoundaryPoint,
183  const PackedBoolList& isTet,
184  const scalarField& cVals,
185  const scalarField& pVals,
186  DynamicList<point>& triPoints,
187  labelList& snappedPoint
188  ) const;
189 
190  //- Generate single point by interpolation or snapping
191  template<class Type>
192  Type generatePoint
193  (
194  const DynamicList<Type>& snappedPoints,
195  const scalar s0,
196  const Type& p0,
197  const label p0Index,
198  const scalar s1,
199  const Type& p1,
200  const label p1Index
201  ) const;
202 
203  template<class Type>
204  void generateTriPoints
205  (
206  const DynamicList<Type>& snapped,
207  const scalar s0,
208  const Type& p0,
209  const label p0Index,
210  const scalar s1,
211  const Type& p1,
212  const label p1Index,
213  const scalar s2,
214  const Type& p2,
215  const label p2Index,
216  const scalar s3,
217  const Type& p3,
218  const label p3Index,
220  ) const;
221 
222  template<class Type>
223  void generateTriPoints
224  (
225  const scalarField& cVals,
226  const scalarField& pVals,
227 
228  const Field<Type>& cCoords,
229  const Field<Type>& pCoords,
230 
231  const DynamicList<Type>& snappedPoints,
232  const labelList& snappedCc,
233  const labelList& snappedPoint,
234 
235  DynamicList<Type>& triPoints,
236  DynamicList<label>& triMeshCells
237  ) const;
238 
239  triSurface stitchTriPoints
240  (
241  const bool checkDuplicates,
242  const List<point>& triPoints,
243  labelList& triPointReverseMap, // unmerged to merged point
244  labelList& triMap // merged to unmerged triangle
245  ) const;
246 
247  //- Check single triangle for (topological) validity
248  static bool validTri(const triSurface&, const label);
249 
250  //- Determine edge-face addressing
251  void calcAddressing
252  (
253  const triSurface& surf,
255  labelList& edgeFace0,
256  labelList& edgeFace1,
257  Map<labelList>& edgeFacesRest
258  ) const;
259 
260  //- Determine orientation
261  static void walkOrientation
262  (
263  const triSurface& surf,
265  const labelList& edgeFace0,
266  const labelList& edgeFace1,
267  const label seedTriI,
268  labelList& flipState
269  );
270 
271  //- Orient surface
272  static void orientSurface
273  (
274  triSurface&,
276  const labelList& edgeFace0,
277  const labelList& edgeFace1,
278  const Map<labelList>& edgeFacesRest
279  );
280 
281  //- Is triangle (given by 3 edges) not fully connected?
282  static bool danglingTriangle
283  (
284  const FixedList<label, 3>& fEdges,
285  const labelList& edgeFace1
286  );
287 
288  //- Mark all non-fully connected triangles
289  static label markDanglingTriangles
290  (
292  const labelList& edgeFace0,
293  const labelList& edgeFace1,
294  const Map<labelList>& edgeFacesRest,
295  boolList& keepTriangles
296  );
297 
298  static triSurface subsetMesh
299  (
300  const triSurface& s,
301  const labelList& newToOldFaces,
302  labelList& oldToNewPoints,
303  labelList& newToOldPoints
304  );
305 
306  //- Combine all triangles inside a cell into a minimal triangulation
307  void combineCellTriangles();
308 
309 public:
310 
311  //- Runtime type information
312  TypeName("isoSurfaceCell");
313 
314 
315  // Constructors
316 
317  //- Construct from dictionary
319  (
320  const polyMesh& mesh,
321  const scalarField& cellValues,
322  const scalarField& pointValues,
323  const scalar iso,
324  const bool regularise,
325  const scalar mergeTol = 1E-6 // fraction of bounding box
326  );
327 
328 
329  // Member Functions
330 
331  //- For every face original cell in mesh
332  const labelList& meshCells() const
333  {
334  return meshCells_;
335  }
336 
337  //- For every unmerged triangle point the point in the triSurface
339  {
340  return triPointMergeMap_;
341  }
342 
343 
344  //- Interpolates cCoords,pCoords. Takes the original fields
345  // used to create the iso surface.
346  template <class Type>
348  (
349  const scalarField& cVals,
350  const scalarField& pVals,
351  const Field<Type>& cCoords,
352  const Field<Type>& pCoords
353  ) const;
354 };
355 
356 
357 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
358 
359 } // End namespace Foam
360 
361 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
362 
363 #ifdef NoRepository
364 # include "isoSurfaceCellTemplates.C"
365 #endif
366 
367 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
368 
369 #endif
370 
371 // ************************ vim: set sw=4 sts=4 et: ************************ //