1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
29 Description of cuts across cells.
31 Description of cut is given as list of vertices and
32 list of edges to be cut (and position on edge).
34 Does some checking of correctness/non-overlapping of cuts. 2x2x2
35 refinement has to be done in three passes since cuts can not overlap
36 (would make addressing too complicated)
38 Introduces concept of 'cut' which is either an existing vertex
42 -# list of cut vertices and list of cut edges. Constructs cell
43 circumference walks ('cellLoops').
45 -# list of cell circumference walks. Will filter them so they
48 -# cellWalker and list of cells to refine (refineCell). Constructs
49 cellLoops and does B. cellWalker is class which can cut a single
50 cell using a plane through the cell centre and in a certain normal
53 CellCuts constructed from cellLoops (B, C) can have multiple cut-edges
54 and/or cut-point as long as there is per face only one (or none) cut across
55 a face, i.e. a face can only be split into two faces.
57 The information available after construction:
58 - pointIsCut, edgeIsCut.
59 - faceSplitCut : the cross-cut of a face.
60 - cellLoops : the loop which cuts across a cell
61 - cellAnchorPoints : per cell the vertices on one side which are
62 considered the anchor points.
64 AnchorPoints: connected loops have to be oriented in the same way to
65 be able to grow new internal faces out of the same bottom faces.
66 (limitation of the mapping procedure). The loop is cellLoop is oriented
67 such that the normal of it points towards the anchorPoints.
68 This is the only routine which uses geometric information.
72 - cut description is very precise. Hard to get right.
73 - do anchor points need to be unique per cell? Very inefficient.
74 - is orientation guaranteed to be correct? Uses geometric info so can go
75 wrong on highly distorted meshes.
76 - is memory inefficient. Probably need to use Maps instead of
82 \*---------------------------------------------------------------------------*/
87 #include "edgeVertex.H"
88 #include "labelList.H"
90 #include "scalarField.H"
91 #include "pointField.H"
92 #include "DynamicList.H"
95 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
100 // Forward declaration of classes
106 /*---------------------------------------------------------------------------*\
107 Class cellCuts Declaration
108 \*---------------------------------------------------------------------------*/
116 // Per point/edge status
118 //- Is mesh point cut
119 boolList pointIsCut_;
124 //- If edge is cut gives weight (0->start() to 1->end())
125 scalarField edgeWeight_;
130 //- Cuts per existing face (includes those along edge of face)
131 // Cuts in no particular order.
132 mutable labelListList* faceCutsPtr_;
134 //- Per face : cut across edge (so not along existing edge)
135 // (can only be one per face)
136 Map<edge> faceSplitCut_;
139 // Cell-cut addressing
141 //- Loop across cell circumference
142 labelListList cellLoops_;
144 //- Number of valid loops in cellLoops_
147 //- For each cut cell the points on the 'anchor' side of the cell.
148 labelListList cellAnchorPoints_;
151 // Private Static Functions
153 //- Find value in first n elements of list.
154 static label findPartIndex
161 //- Create boolList with all labels specified set to true
162 // (and rest to false)
163 static boolList expand(const label size, const labelList& labels);
165 //- Create scalarField with all specified labels set to corresponding
166 // value in scalarField.
167 static scalarField expand
174 //- Returns -1 or index of first element of lst that cannot be found
176 static label firstUnique
178 const labelList& lst,
182 // Private Member Functions
184 //- Debugging: write cell's edges and any cut vertices and edges
185 // (so no cell loop determined yet)
186 void writeUncutOBJ(const fileName&, const label cellI) const;
188 //- Debugging: write cell's edges, loop and anchors to directory.
193 const pointField& loopPoints,
194 const labelList& anchors
197 //- Find face on cell using the two edges.
206 //- Find face on cell using an edge and a vertex.
207 label edgeVertexToFace
214 //- Find face using two vertices (guaranteed not to be along edge)
215 label vertexVertexToFace
225 //- Calculate faceCuts in face vertex order.
226 void calcFaceCuts() const;
229 // Loop (cuts on cell circumference) calculation
231 //- Find edge (or -1) on faceI using vertices v0,v1
239 //- Find face on which all cuts are (very rare) or -1.
240 label loopFace(const label cellI, const labelList& loop) const;
242 //- Cross otherCut into next faces (not exclude0, exclude1)
246 const label startCut,
248 const label exclude0,
249 const label exclude1,
251 const label otherCut,
257 //- Cross cut (which is edge on faceI) onto next face
261 const label startCut,
263 const label otherCut,
269 // wrapper around visited[nVisited++] = cut. Checks for duplicate
279 //- Walk across faceI following cuts, starting at cut. Stores cuts
284 const label startCut,
289 label& beforeLastCut,
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
300 const label startCut, // overall starting cut
302 const label prevCut, // current cut
307 //- Determine for every cut cell the face it is cut by.
308 void calcCellLoops(const labelList& cutCells);
313 //- Are there enough faces on anchor side of cellI?
317 const labelList& anchorPoints
320 //- Walk unset edges of single cell from starting point and
321 // marks visited edges and vertices with status.
328 Map<label>& edgeStatus,
329 Map<label>& pointStatus
332 //- Check anchor points on 'outside' of loop
333 bool loopAnchorConsistent
336 const pointField& loopPts,
337 const labelList& anchorPoints
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.
346 const labelList& loop,
347 const pointField& loopPts,
349 labelList& anchorPoints
352 //- Returns coordinates of points on loop with explicitly provided
354 pointField loopPoints
356 const labelList& loop,
357 const scalarField& loopWeights
360 //- Returns weights of loop. Inverse of loopPoints.
361 scalarField loopWeights(const labelList& loop) const;
363 //- Check if cut edges in loop are compatible with ones in
367 const labelList& loop,
368 const scalarField& loopWeights
371 //- Counts number of cuts on face.
375 const labelList& loop
378 //- Determines if loop through cellI consistent with existing
380 bool conservativeValidLoop
383 const labelList& loop
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.
393 const labelList& loop,
394 const scalarField& loopWeights,
395 Map<edge>& newFaceSplitCut,
396 labelList& anchorPoints
399 //- Update basic cut information from cellLoops. Assumes cellLoops_
400 // already set and consistent.
401 void setFromCellLoops();
403 //- Update basic cut information for single cell from cellLoop.
407 const labelList& loop,
408 const scalarField& loopWeights
411 //- Update basic cut information from cellLoops. Checks for
412 // consistency with existing cut pattern.
413 void setFromCellLoops
415 const labelList& cellLabels,
416 const labelListList& cellLoops,
417 const List<scalarField>& cellLoopWeights
420 //- Cut cells and update basic cut information from cellLoops.
421 // Checks each loop for consistency with existing cut pattern.
422 void setFromCellCutter
425 const List<refineCell>& refCells
428 //- Same as above but now cut with prescribed plane.
429 void setFromCellCutter
432 const labelList& cellLabels,
436 //- Set orientation of loops
437 void orientPlanesAndLoops();
439 //- top level driver: adressing calculation and loop detection
440 void calcLoopsAndAddressing(const labelList& cutCells);
442 //- Check various consistencies.
446 //- Disallow default bitwise copy construct
447 cellCuts(const cellCuts&);
449 //- Disallow default bitwise assignment
450 void operator=(const cellCuts&);
455 //- Runtime type information
456 ClassName("cellCuts");
461 //- Construct from cells to cut and pattern of cuts
464 const polyMesh& mesh,
465 const labelList& cutCells,
466 const labelList& meshVerts,
467 const labelList& meshEdges,
468 const scalarField& meshEdgeWeights
471 //- Construct from pattern of cuts. Detect cells to cut.
474 const polyMesh& mesh,
475 const labelList& meshVerts,
476 const labelList& meshEdges,
477 const scalarField& meshEdgeWeights
480 //- Construct from complete cellLoops through specified cells.
481 // Checks for consistency.
482 // Constructs cut-cut addressing and cellAnchorPoints.
485 const polyMesh& mesh,
486 const labelList& cellLabels,
487 const labelListList& cellLoops,
488 const List<scalarField>& cellEdgeWeights
491 //- Construct from list of cells to cut and direction to cut in
492 // (always through cell centre) and celllooper.
495 const polyMesh& mesh,
496 const cellLooper& cellCutter,
497 const List<refineCell>& refCells
500 //- Construct from list of cells to cut and plane to cut with and
501 // celllooper. (constructor above always cuts through cell centre)
504 const polyMesh& mesh,
505 const cellLooper& cellCutter,
506 const labelList& cellLabels,
507 const List<plane>& cutPlanes
510 //- Construct from components
513 const polyMesh& mesh,
514 const boolList& pointIsCut,
515 const boolList& edgeIsCut,
516 const scalarField& edgeWeight,
517 const Map<edge>& faceSplitCut,
518 const labelListList& cellLoops,
520 const labelListList& cellAnchorPoints
528 //- Clear out demand driven storage
536 //- Is mesh point cut
537 const boolList& pointIsCut() const
543 const boolList& edgeIsCut() const
548 //- If edge is cut gives weight (ratio between start() and end())
549 const scalarField& edgeWeight() const
554 //- Cuts per existing face (includes those along edge of face)
555 // Cuts in no particular order
556 const labelListList& faceCuts() const
562 return *faceCutsPtr_;
565 //- Gives for split face the two cuts that split the face into two.
566 const Map<edge>& faceSplitCut() const
568 return faceSplitCut_;
571 //- For each cut cell the cut along the circumference.
572 const labelListList& cellLoops() const
577 //- Number of valid cell loops
583 //- For each cut cell the points on the 'anchor' side of the cell.
584 const labelListList& cellAnchorPoints() const
586 return cellAnchorPoints_;
592 //- Returns coordinates of points on loop for given cell.
593 // Uses cellLoops_ and edgeWeight_
594 pointField loopPoints(const label cellI) const;
596 //- Invert anchor point selection.
597 labelList nonAnchorPoints
599 const labelList& cellPoints,
600 const labelList& anchorPoints,
601 const labelList& loop
604 //- Flip loop for cellI. Updates anchor points as well.
605 void flip(const label cellI);
607 //- Flip loop for cellI. Does not update anchors. Use with care
608 // (only if you're sure loop orientation is wrong)
609 void flipLoopOnly(const label cellI);
614 //- debugging:Write list of cuts to stream in OBJ format
618 const pointField& loopPoints,
622 //- debugging:Write all of cuts to stream in OBJ format
623 void writeOBJ(Ostream& os) const;
625 //- debugging:Write edges of cell and loop
626 void writeCellOBJ(const fileName& dir, const label cellI) const;
631 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
633 } // End namespace Foam
635 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
639 // ************************************************************************* //