cellCuts.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::cellCuts
26 
27 Description
28  Description of cuts across cells.
29 
30  Description of cut is given as list of vertices and
31  list of edges to be cut (and position on edge).
32 
33  Does some checking of correctness/non-overlapping of cuts. 2x2x2
34  refinement has to be done in three passes since cuts can not overlap
35  (would make addressing too complicated)
36 
37  Introduces concept of 'cut' which is either an existing vertex
38  or a edge.
39 
40  Input can either be
41  -# list of cut vertices and list of cut edges. Constructs cell
42  circumference walks ('cellLoops').
43 
44  -# list of cell circumference walks. Will filter them so they
45  don't overlap.
46 
47  -# cellWalker and list of cells to refine (refineCell). Constructs
48  cellLoops and does B. cellWalker is class which can cut a single
49  cell using a plane through the cell centre and in a certain normal
50  direction
51 
52  CellCuts constructed from cellLoops (B, C) can have multiple cut-edges
53  and/or cut-point as long as there is per face only one (or none) cut across
54  a face, i.e. a face can only be split into two faces.
55 
56  The information available after construction:
57  - pointIsCut, edgeIsCut.
58  - faceSplitCut : the cross-cut of a face.
59  - cellLoops : the loop which cuts across a cell
60  - cellAnchorPoints : per cell the vertices on one side which are
61  considered the anchor points.
62 
63  AnchorPoints: connected loops have to be oriented in the same way to
64  be able to grow new internal faces out of the same bottom faces.
65  (limitation of the mapping procedure). The loop is cellLoop is oriented
66  such that the normal of it points towards the anchorPoints.
67  This is the only routine which uses geometric information.
68 
69 
70  Limitations:
71  - cut description is very precise. Hard to get right.
72  - do anchor points need to be unique per cell? Very inefficient.
73  - is orientation guaranteed to be correct? Uses geometric info so can go
74  wrong on highly distorted meshes.
75  - is memory inefficient. Probably need to use Maps instead of
76  labelLists.
77 
78 SourceFiles
79  cellCuts.C
80 
81 \*---------------------------------------------------------------------------*/
82 
83 #ifndef cellCuts_H
84 #define cellCuts_H
85 
86 #include "edgeVertex.H"
87 #include "labelList.H"
88 #include "boolList.H"
89 #include "scalarField.H"
90 #include "pointField.H"
91 #include "DynamicList.H"
92 #include "typeInfo.H"
93 
94 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
95 
96 namespace Foam
97 {
98 
99 // Forward declaration of classes
100 class polyMesh;
101 class cellLooper;
102 class refineCell;
103 class plane;
104 
105 /*---------------------------------------------------------------------------*\
106  Class cellCuts Declaration
107 \*---------------------------------------------------------------------------*/
109 class cellCuts
110 :
111  public edgeVertex
112 {
113  // Private data
114 
115  // Per point/edge status
116 
117  //- Is mesh point cut
118  boolList pointIsCut_;
119 
120  //- Is edge cut
121  boolList edgeIsCut_;
122 
123  //- If edge is cut gives weight (0->start() to 1->end())
124  scalarField edgeWeight_;
125 
126 
127  // Cut addressing
128 
129  //- Cuts per existing face (includes those along edge of face)
130  // Cuts in no particular order.
131  mutable labelListList* faceCutsPtr_;
132 
133  //- Per face : cut across edge (so not along existing edge)
134  // (can only be one per face)
135  Map<edge> faceSplitCut_;
136 
137 
138  // Cell-cut addressing
139 
140  //- Loop across cell circumference
141  labelListList cellLoops_;
142 
143  //- Number of valid loops in cellLoops_
144  label nLoops_;
145 
146  //- For each cut cell the points on the 'anchor' side of the cell.
147  labelListList cellAnchorPoints_;
148 
149 
150  // Private Static Functions
151 
152  //- Find value in first n elements of list.
153  static label findPartIndex
154  (
155  const labelList&,
156  const label n,
157  const label val
158  );
159 
160  //- Create boolList with all labels specified set to true
161  // (and rest to false)
162  static boolList expand(const label size, const labelList& labels);
163 
164  //- Create scalarField with all specified labels set to corresponding
165  // value in scalarField.
166  static scalarField expand
167  (
168  const label,
169  const labelList&,
170  const scalarField&
171  );
172 
173  //- Returns -1 or index of first element of lst that cannot be found
174  // in map.
175  static label firstUnique
176  (
177  const labelList& lst,
178  const Map<label>&
179  );
180 
181  // Private Member Functions
182 
183  //- Debugging: write cell's edges and any cut vertices and edges
184  // (so no cell loop determined yet)
185  void writeUncutOBJ(const fileName&, const label celli) const;
186 
187  //- Debugging: write cell's edges, loop and anchors to directory.
188  void writeOBJ
189  (
190  const fileName& dir,
191  const label celli,
192  const pointField& loopPoints,
193  const labelList& anchors
194  ) const;
195 
196  //- Find face on cell using the two edges.
197  label edgeEdgeToFace
198  (
199  const label celli,
200  const label edgeA,
201  const label edgeB
202  ) const;
203 
204 
205  //- Find face on cell using an edge and a vertex.
206  label edgeVertexToFace
207  (
208  const label celli,
209  const label edgeI,
210  const label vertI
211  ) const;
212 
213  //- Find face using two vertices (guaranteed not to be along edge)
214  label vertexVertexToFace
215  (
216  const label celli,
217  const label vertA,
218  const label vertB
219  ) const;
220 
221 
222  // Cut addressing
223 
224  //- Calculate faceCuts in face vertex order.
225  void calcFaceCuts() const;
226 
227 
228  // Loop (cuts on cell circumference) calculation
229 
230  //- Find edge (or -1) on facei using vertices v0,v1
231  label findEdge
232  (
233  const label facei,
234  const label v0,
235  const label v1
236  ) const;
237 
238  //- Find face on which all cuts are (very rare) or -1.
239  label loopFace(const label celli, const labelList& loop) const;
240 
241  //- Cross otherCut into next faces (not exclude0, exclude1)
242  bool walkPoint
243  (
244  const label celli,
245  const label startCut,
246 
247  const label exclude0,
248  const label exclude1,
249 
250  const label otherCut,
251 
252  label& nVisited,
253  labelList& visited
254  ) const;
255 
256  //- Cross cut (which is edge on facei) onto next face
257  bool crossEdge
258  (
259  const label celli,
260  const label startCut,
261  const label facei,
262  const label otherCut,
263 
264  label& nVisited,
265  labelList& visited
266  ) const;
267 
268  // wrapper around visited[nVisited++] = cut. Checks for duplicate
269  // cuts.
270  bool addCut
271  (
272  const label celli,
273  const label cut,
274  label& nVisited,
275  labelList& visited
276  ) const;
277 
278  //- Walk across facei following cuts, starting at cut. Stores cuts
279  // visited
280  // Returns true if valid walk.
281  bool walkFace
282  (
283  const label celli,
284  const label startCut,
285  const label facei,
286  const label cut,
287 
288  label& lastCut,
289  label& beforeLastCut,
290  label& nVisited,
291  labelList& visited
292  ) const;
293 
294  //- Walk across cuts (cut edges or cut vertices) of cell. Stops when
295  // hit cut already visited. Returns true when loop of 3 or more
296  // vertices found.
297  bool walkCell
298  (
299  const label celli,
300  const label startCut, // overall starting cut
301  const label facei,
302  const label prevCut, // current cut
303  label& nVisited,
304  labelList& visited
305  ) const;
306 
307  //- Determine for every cut cell the face it is cut by.
308  void calcCellLoops(const labelList& cutCells);
309 
310 
311  // Cell anchoring
312 
313  //- Are there enough faces on anchor side of celli?
314  bool checkFaces
315  (
316  const label celli,
317  const labelList& anchorPoints
318  ) const;
319 
320  //- Walk unset edges of single cell from starting point and
321  // marks visited edges and vertices with status.
322  void walkEdges
323  (
324  const label celli,
325  const label pointi,
326  const label status,
327 
328  Map<label>& edgeStatus,
329  Map<label>& pointStatus
330  ) const;
331 
332  //- Check anchor points on 'outside' of loop
333  bool loopAnchorConsistent
334  (
335  const label celli,
336  const pointField& loopPts,
337  const labelList& anchorPoints
338  ) const;
339 
340  //- Determines set of anchor points given a loop. The loop should
341  // split the cell into two. Returns true if a valid set of anchor
342  // points determined, false otherwise.
343  bool calcAnchors
344  (
345  const label celli,
346  const labelList& loop,
347  const pointField& loopPts,
348 
349  labelList& anchorPoints
350  ) const;
351 
352  //- Returns coordinates of points on loop with explicitly provided
353  // weights.
354  pointField loopPoints
355  (
356  const labelList& loop,
357  const scalarField& loopWeights
358  ) const;
359 
360  //- Returns weights of loop. Inverse of loopPoints.
361  scalarField loopWeights(const labelList& loop) const;
362 
363  //- Check if cut edges in loop are compatible with ones in
364  // edgeIsCut_
365  bool validEdgeLoop
366  (
367  const labelList& loop,
368  const scalarField& loopWeights
369  ) const;
370 
371  //- Counts number of cuts on face.
372  label countFaceCuts
373  (
374  const label facei,
375  const labelList& loop
376  ) const;
377 
378  //- Determine compatibility of loop with existing cut pattern.
379  // Does not use cut-addressing (faceCuts_, cutCuts_)
380  bool conservativeValidLoop
381  (
382  const label celli,
383  const labelList& loop
384  ) const;
385 
386  //- Check if loop is compatible with existing cut pattern in
387  // pointIsCut, edgeIsCut, faceSplitCut.
388  // Calculates and returns for current cell the cut faces and the
389  // points on one side of the loop.
390  bool validLoop
391  (
392  const label celli,
393  const labelList& loop,
394  const scalarField& loopWeights,
395  Map<edge>& newFaceSplitCut,
396  labelList& anchorPoints
397  ) const;
398 
399  //- Update basic cut information from cellLoops.
400  // Assumes cellLoops_ and edgeWeight_ already set and consistent.
401  // Does not use any other information.
402  void setFromCellLoops();
403 
404  //- Update basic cut information for single cell from cellLoop.
405  bool setFromCellLoop
406  (
407  const label celli,
408  const labelList& loop,
409  const scalarField& loopWeights
410  );
411 
412  //- Update basic cut information from cellLoops. Checks for
413  // consistency with existing cut pattern.
414  void setFromCellLoops
415  (
416  const labelList& cellLabels,
417  const labelListList& cellLoops,
418  const List<scalarField>& cellLoopWeights
419  );
420 
421  //- Cut cells and update basic cut information from cellLoops.
422  // Checks each loop for consistency with existing cut pattern.
423  void setFromCellCutter
424  (
425  const cellLooper&,
426  const List<refineCell>& refCells
427  );
428 
429  //- Same as above but now cut with prescribed plane (instead of
430  // just normal).
431  void setFromCellCutter
432  (
433  const cellLooper&,
434  const labelList& cellLabels,
435  const List<plane>&
436  );
437 
438  //- Set orientation of loops
439  void orientPlanesAndLoops();
440 
441  //- Top level driver: adressing calculation and loop detection
442  // (loops splitting cells).
443  void calcLoopsAndAddressing(const labelList& cutCells);
444 
445  //- Check various consistencies.
446  void check() const;
447 
448 
449  //- Disallow default bitwise copy construct
450  cellCuts(const cellCuts&);
451 
452  //- Disallow default bitwise assignment
453  void operator=(const cellCuts&);
454 
455 
456 public:
457 
458  //- Runtime type information
459  ClassName("cellCuts");
460 
461 
462  // Constructors
463 
464  //- Construct from cells to cut and pattern of cuts
465  cellCuts
466  (
467  const polyMesh& mesh,
468  const labelList& cutCells,
469  const labelList& meshVerts,
470  const labelList& meshEdges,
471  const scalarField& meshEdgeWeights
472  );
473 
474  //- Construct from pattern of cuts. Detect cells to cut.
475  cellCuts
476  (
477  const polyMesh& mesh,
478  const labelList& meshVerts,
479  const labelList& meshEdges,
480  const scalarField& meshEdgeWeights
481  );
482 
483  //- Construct from complete cellLoops through specified cells.
484  // Checks for consistency.
485  // Constructs cut-cut addressing and cellAnchorPoints.
486  cellCuts
487  (
488  const polyMesh& mesh,
489  const labelList& cellLabels,
490  const labelListList& cellLoops,
491  const List<scalarField>& cellEdgeWeights
492  );
493 
494  //- Construct from list of cells to cut and direction to cut in
495  // (always through cell centre) and celllooper.
496  cellCuts
497  (
498  const polyMesh& mesh,
499  const cellLooper& cellCutter,
500  const List<refineCell>& refCells
501  );
502 
503  //- Construct from list of cells to cut and plane to cut with and
504  // celllooper. (constructor above always cuts through cell centre)
505  cellCuts
506  (
507  const polyMesh& mesh,
508  const cellLooper& cellCutter,
509  const labelList& cellLabels,
510  const List<plane>& cutPlanes
511  );
512 
513  //- Construct from components
514  cellCuts
515  (
516  const polyMesh& mesh,
517  const boolList& pointIsCut,
518  const boolList& edgeIsCut,
519  const scalarField& edgeWeight,
520  const Map<edge>& faceSplitCut,
521  const labelListList& cellLoops,
522  const label nLoops,
524  );
525 
526 
527  //- Destructor
528  ~cellCuts();
529 
530  //- Clear out demand driven storage
531  void clearOut();
532 
533 
534  // Member Functions
535 
536  // Access
537 
538  //- Is mesh point cut
539  const boolList& pointIsCut() const
540  {
541  return pointIsCut_;
542  }
543 
544  //- Is edge cut
545  const boolList& edgeIsCut() const
546  {
547  return edgeIsCut_;
548  }
549 
550  //- If edge is cut gives weight (ratio between start() and end())
551  const scalarField& edgeWeight() const
552  {
553  return edgeWeight_;
554  }
555 
556  //- Cuts per existing face (includes those along edge of face)
557  // Cuts in no particular order
558  const labelListList& faceCuts() const
559  {
560  if (!faceCutsPtr_)
561  {
562  calcFaceCuts();
563  }
564  return *faceCutsPtr_;
565  }
566 
567  //- Gives for split face the two cuts that split the face into two.
568  const Map<edge>& faceSplitCut() const
569  {
570  return faceSplitCut_;
571  }
572 
573  //- For each cut cell the cut along the circumference.
574  const labelListList& cellLoops() const
575  {
576  return cellLoops_;
577  }
578 
579  //- Number of valid cell loops
580  label nLoops() const
581  {
582  return nLoops_;
583  }
584 
585  //- For each cut cell the points on the 'anchor' side of the cell.
586  const labelListList& cellAnchorPoints() const
587  {
588  return cellAnchorPoints_;
589  }
590 
591 
592  // Other
593 
594  //- Returns coordinates of points on loop for given cell.
595  // Uses cellLoops_ and edgeWeight_
596  pointField loopPoints(const label celli) const;
597 
598  //- Invert anchor point selection.
600  (
601  const labelList& cellPoints,
602  const labelList& anchorPoints,
603  const labelList& loop
604  ) const;
605 
606  //- Flip loop for celli. Updates anchor points as well.
607  void flip(const label celli);
608 
609  //- Flip loop for celli. Does not update anchors. Use with care
610  // (only if you're sure loop orientation is wrong)
611  void flipLoopOnly(const label celli);
612 
613 
614  // Write
615 
616  //- debugging:Write list of cuts to stream in OBJ format
617  void writeOBJ
618  (
619  Ostream& os,
620  const pointField& loopPoints,
621  label& vertI
622  ) const;
623 
624  //- debugging:Write all of cuts to stream in OBJ format
625  void writeOBJ(Ostream& os) const;
626 
627  //- debugging:Write edges of cell and loop
628  void writeCellOBJ(const fileName& dir, const label celli) const;
629 
630 };
631 
632 
633 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
634 
635 } // End namespace Foam
636 
637 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
638 
639 #endif
640 
641 // ************************************************************************* //
const boolList & edgeIsCut() const
Is edge cut.
Definition: cellCuts.H:544
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
A class for handling file names.
Definition: fileName.H:69
ClassName("cellCuts")
Runtime type information.
void flip(const label celli)
Flip loop for celli. Updates anchor points as well.
Definition: cellCuts.C:2958
Abstract base class. Concrete implementations know how to cut a cell (i.e. determine a loop around th...
Definition: cellLooper.H:69
label nLoops() const
Number of valid cell loops.
Definition: cellCuts.H:579
Description of cuts across cells.
Definition: cellCuts.H:108
const polyMesh & mesh() const
Definition: edgeVertex.H:98
const scalarField & edgeWeight() const
If edge is cut gives weight (ratio between start() and end())
Definition: cellCuts.H:550
void clearOut()
Clear out demand driven storage.
Definition: cellCuts.C:2925
Combines edge or vertex in single label. Used to specify cuts across cell circumference.
Definition: edgeVertex.H:52
const boolList & pointIsCut() const
Is mesh point cut.
Definition: cellCuts.H:538
~cellCuts()
Destructor.
Definition: cellCuts.C:2919
const labelListList & cellLoops() const
For each cut cell the cut along the circumference.
Definition: cellCuts.H:573
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const labelListList & faceCuts() const
Cuts per existing face (includes those along edge of face)
Definition: cellCuts.H:557
const Map< edge > & faceSplitCut() const
Gives for split face the two cuts that split the face into two.
Definition: cellCuts.H:567
const labelListList & cellAnchorPoints() const
For each cut cell the points on the &#39;anchor&#39; side of the cell.
Definition: cellCuts.H:585
label n
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
void writeCellOBJ(const fileName &dir, const label celli) const
debugging:Write edges of cell and loop
Definition: cellCuts.C:3021
labelList nonAnchorPoints(const labelList &cellPoints, const labelList &anchorPoints, const labelList &loop) const
Invert anchor point selection.
Definition: cellCuts.C:1167
Namespace for OpenFOAM.
void flipLoopOnly(const label celli)
Flip loop for celli. Does not update anchors. Use with care.
Definition: cellCuts.C:2975
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49