OpenFOAM logo
Open Source CFD Toolkit

cellCuts.H

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2005 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software; you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM; if not, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     cellCuts
00027 
00028 Description
00029     Description of cuts across cells.
00030 
00031     Description of cut is given as list of vertices and
00032     list of edges to be cut (and position on edge).
00033 
00034     Does some checking of correctness/non-overlapping of cuts. 2x2x2
00035     refinement has to be done in three passes since cuts can not overlap
00036     (would make addressing too complicated)
00037 
00038     Introduces concept of 'cut' which is either an existing vertex
00039     or a edge.
00040 
00041     Input can either be
00042     A. list of cut vertices and list of cut edges. Constructs cell circumference
00043     walks ('cellLoops').
00044 
00045     B. list of cell circumference walks. Will filter them so they
00046     don't overlap.
00047 
00048     C. cellWalker and list of cells to refine (refineCell). Constructs cellLoops
00049     and does B. cellWalker is class which can cut a single cell using a
00050     plane through the cell centre and in a certain normal direction
00051 
00052 
00053     CellCuts constructed from cellLoops (B, C) can have multiple cut-edges
00054     and/or cut-point as long as there is per face only one (or none) cut across
00055     a face, i.e. a face can only be split into two faces.
00056 
00057     The information available after construction:
00058     - pointIsCut, edgeIsCut.
00059     - faceSplitCut : the cross-cut of a face.
00060     - cellLoops : the loop which cuts across a cell
00061     - cellAnchorPoints : per cell the vertices on one side which are considered
00062         the anchor points.
00063 
00064     AnchorPoints: connected loops have to be oriented in the same way to 
00065     be able to grow new internal faces out of the same bottom faces.
00066     (limitation of the mapping procedure). The loop is cellLoop is oriented
00067     such that the normal of it points towards the anchorPoints.
00068     This is the only routine which uses geometric information.
00069 
00070 
00071     Limitations:
00072     - cut description is very precise. Hard to get right.
00073     - do anchor points need to be unique per cell? Very inefficient.
00074     - is orientation guaranteed to be correct? Uses geometric info so can go
00075       wrong on highly distorted meshes.
00076     - is memory inefficient. Probably need to use Maps instead of
00077       labelLists.
00078 
00079 SourceFiles
00080     cellCuts.C
00081 
00082 \*---------------------------------------------------------------------------*/
00083 
00084 #ifndef cellCuts_H
00085 #define cellCuts_H
00086 
00087 #include "edgeVertex.H"
00088 #include "labelList.H"
00089 #include "boolList.H"
00090 #include "scalarField.H"
00091 #include "pointField.H"
00092 #include "DynamicList.H"
00093 #include "typeInfo.H"
00094 
00095 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00096 
00097 namespace Foam
00098 {
00099 
00100 // Class forward declarations
00101 class polyMesh;
00102 class cellLooper;
00103 class refineCell;
00104 class plane;
00105 
00106 /*---------------------------------------------------------------------------*\
00107                            Class cellCuts Declaration
00108 \*---------------------------------------------------------------------------*/
00109 
00110 class cellCuts
00111 :
00112     public edgeVertex
00113 {
00114     // Private data
00115 
00116         // Per point/edge status
00117 
00118             //- Is mesh point cut
00119             boolList pointIsCut_;
00120 
00121             //- Is edge cut
00122             boolList edgeIsCut_;
00123 
00124             //- If edge is cut gives weight (0->start() to 1->end())
00125             scalarField edgeWeight_;
00126 
00127 
00128         // Cut addressing
00129 
00130             //- Cuts per existing face (includes those along edge of face)
00131             //  Cuts in no particular order.
00132             mutable labelListList* faceCutsPtr_;
00133 
00134             //- Per face : cut across edge (so not along existing edge)
00135             //  (can only be one per face)
00136             Map<edge> faceSplitCut_;
00137 
00138 
00139         // Cell-cut addressing
00140 
00141             //- Loop across cell circumference
00142             labelListList cellLoops_;
00143 
00144             //- Number of valid loops in cellLoops_
00145             label nLoops_;
00146 
00147             //- For each cut cell the points on the 'anchor' side of the cell.
00148             labelListList cellAnchorPoints_;
00149 
00150 
00151     // Private Static Functions
00152 
00153         //- Find value in first n elements of list.
00154         static label findPartIndex
00155         (
00156             const labelList&,
00157             const label n,
00158             const label val
00159         );
00160 
00161         //- Create boolList with all labels specified set to true
00162         //  (and rest to false)
00163         static boolList expand(const label size, const labelList& labels);
00164 
00165         //- Create scalarField with all specified labels set to corresponding
00166         //  value in scalarField.
00167         static scalarField expand
00168         (
00169             const label,
00170             const labelList&,
00171             const scalarField&
00172         );
00173 
00174         //- Returns -1 or index of first element of lst that cannot be found
00175         //  in map.
00176         static label firstUnique
00177         (
00178             const labelList& lst,
00179             const Map<label>&
00180         );
00181 
00182     // Private Member Functions
00183 
00184         //- Debugging: write cell's edges and any cut vertices and edges
00185         //  (so no cell loop determined yet)
00186         void writeUncutOBJ(const fileName&, const label cellI) const;
00187 
00188         //- Debugging: write cell's edges, loop and anchors to directory.
00189         void writeOBJ
00190         (
00191             const fileName& dir,
00192             const label cellI,
00193             const pointField& loopPoints,
00194             const labelList& anchors
00195         ) const;
00196 
00197         //- Find face on cell using the two edges.
00198         label edgeEdgeToFace
00199         (
00200             const label cellI,
00201             const label edgeA,
00202             const label edgeB
00203         ) const;
00204 
00205 
00206         //- Find face on cell using an edge and a vertex.
00207         label edgeVertexToFace
00208         (
00209             const label cellI,
00210             const label edgeI,
00211             const label vertI
00212         ) const;
00213 
00214         //- Find face using two vertices (guaranteed not to be along edge)
00215         label vertexVertexToFace
00216         (
00217             const label cellI,
00218             const label vertA,
00219             const label vertB
00220         ) const;
00221 
00222 
00223         // Cut addressing
00224 
00225             //- Calculate faceCuts in face vertex order.
00226             void calcFaceCuts() const;
00227 
00228 
00229         // Loop (cuts on cell circumference) calculation
00230 
00231             //- Find edge (or -1) on faceI using vertices v0,v1
00232             label findEdge
00233             (
00234                 const label faceI,
00235                 const label v0,
00236                 const label v1
00237             ) const;
00238 
00239             //- Find face on which all cuts are (very rare) or -1.
00240             label loopFace(const label cellI, const labelList& loop) const;
00241 
00242             //- Cross otherCut into next faces (not exclude0, exclude1)
00243             bool walkPoint
00244             (
00245                 const label cellI,
00246                 const label startCut,
00247 
00248                 const label exclude0,
00249                 const label exclude1,
00250 
00251                 const label otherCut,
00252 
00253                 label& nVisited,
00254                 labelList& visited
00255             ) const;
00256 
00257             //- Cross cut (which is edge on faceI) onto next face
00258             bool crossEdge
00259             (
00260                 const label cellI,
00261                 const label startCut,
00262                 const label faceI,
00263                 const label otherCut,
00264 
00265                 label& nVisited,
00266                 labelList& visited
00267             ) const;
00268 
00269             // wrapper around visited[nVisited++] = cut. Checks for duplicate
00270             // cuts.
00271             bool addCut
00272             (
00273                 const label cellI,
00274                 const label cut,
00275                 label& nVisited,
00276                 labelList& visited
00277             ) const;
00278 
00279             //- Walk across faceI following cuts, starting at cut. Stores cuts
00280             //  visited
00281             bool walkFace
00282             (
00283                 const label cellI,
00284                 const label startCut,
00285                 const label faceI,
00286                 const label cut,
00287 
00288                 label& lastCut,
00289                 label& beforeLastCut,
00290                 label& nVisited,
00291                 labelList& visited
00292             ) const;
00293 
00294             //- Walk across cuts (cut edges or cut vertices) of cell. Stops when
00295             //  hit cut  already visited. Returns true when loop of 3 or more
00296             //  vertices found.
00297             bool walkCell
00298             (
00299                 const label cellI,
00300                 const label startCut,   // overall starting cut
00301                 const label faceI,
00302                 const label prevCut,    // current cut
00303                 label& nVisited,
00304                 labelList& visited
00305             ) const;
00306 
00307             //- Determine for every cut cell the face it is cut by.
00308             void calcCellLoops(const labelList& cutCells);
00309 
00310 
00311         // Cell anchoring
00312 
00313             //- Are there enough faces on anchor side of cellI?
00314             bool checkFaces
00315             (
00316                 const label cellI,
00317                 const labelList& anchorPoints
00318             ) const;
00319 
00320             //- Walk unset edges of single cell from starting point and
00321             //  marks visited edges and vertices with status.
00322             void walkEdges
00323             (
00324                 const label cellI,
00325                 const label pointI,
00326                 const label status,
00327 
00328                 Map<label>& edgeStatus,
00329                 Map<label>& pointStatus
00330             ) const;
00331 
00332             //- Check anchor points on 'outside' of loop
00333             bool loopAnchorConsistent
00334             (
00335                 const label cellI,
00336                 const pointField& loopPts,
00337                 const labelList& anchorPoints
00338             ) const;
00339 
00340             //- Determines set of anchor points given a loop. The loop should
00341             //  split the cell into two. Returns true if a valid set of anchor
00342             //  points determined, false otherwise.
00343             bool calcAnchors
00344             (
00345                 const label cellI,
00346                 const labelList& loop,
00347                 const pointField& loopPts,
00348 
00349                 labelList& anchorPoints
00350             ) const;
00351 
00352             //- Returns coordinates of points on loop with explicitly provided
00353             //  weights.
00354             pointField loopPoints
00355             (
00356                 const labelList& loop,
00357                 const scalarField& loopWeights
00358             ) const;
00359 
00360             //- Returns weights of loop. Inverse of loopPoints.
00361             scalarField loopWeights(const labelList& loop) const;
00362 
00363             //- Check if cut edges in loop are compatible with ones in
00364             //  edgeIsCut_
00365             bool validEdgeLoop
00366             (
00367                 const labelList& loop,
00368                 const scalarField& loopWeights
00369             ) const;
00370 
00371             //- Counts number of cuts on face.
00372             label countFaceCuts
00373             (
00374                 const label faceI,
00375                 const labelList& loop
00376             ) const;
00377 
00378             //- Determines if loop through cellI consistent with existing
00379             //  pattern.
00380             bool conservativeValidLoop
00381             (
00382                 const label cellI,
00383                 const labelList& loop
00384             ) const;
00385 
00386             //- Check if loop is compatible with existing cut pattern in
00387             //  pointIsCut, edgeIsCut, faceSplitCut.
00388             //  Calculates and returns for current cell the cut faces and the
00389             //  points on one side of the loop.
00390             bool validLoop
00391             (
00392                 const label cellI,
00393                 const labelList& loop,
00394                 const scalarField& loopWeights,
00395                 Map<edge>& newFaceSplitCut,
00396                 labelList& anchorPoints
00397             ) const;
00398 
00399             //- Update basic cut information from cellLoops. Assumes cellLoops_
00400             //  already set and consistent.
00401             void setFromCellLoops();
00402 
00403             //- Update basic cut information for single cell from cellLoop.
00404             bool setFromCellLoop
00405             (
00406                 const label cellI,
00407                 const labelList& loop,
00408                 const scalarField& loopWeights
00409             );
00410 
00411             //- Update basic cut information from cellLoops. Checks for
00412             //  consistency with existing cut pattern.
00413             void setFromCellLoops
00414             (
00415                 const labelList& cellLabels,
00416                 const labelListList& cellLoops,
00417                 const List<scalarField>& cellLoopWeights
00418             );
00419 
00420             //- Cut cells and update basic cut information from cellLoops.
00421             //  Checks each loop for consistency with existing cut pattern.
00422             void setFromCellCutter
00423             (
00424                 const cellLooper&,
00425                 const List<refineCell>& refCells
00426             );
00427 
00428             //- Same as above but now cut with prescribed plane.
00429             void setFromCellCutter
00430             (
00431                 const cellLooper&,
00432                 const labelList& cellLabels,
00433                 const List<plane>&
00434             );
00435 
00436             //- Set orientation of loops
00437             void orientPlanesAndLoops();
00438             
00439             //- top level driver: adressing calculation and loop detection
00440             void calcLoopsAndAddressing(const labelList& cutCells);
00441 
00442             //- Check various consistencies.
00443             void check() const;
00444 
00445 
00446         //- Disallow default bitwise copy construct
00447         cellCuts(const cellCuts&);
00448 
00449         //- Disallow default bitwise assignment
00450         void operator=(const cellCuts&);
00451 
00452 
00453 public:
00454 
00455     //- Runtime type information
00456     ClassName("cellCuts");
00457 
00458 
00459     // Constructors
00460 
00461         //- Construct from cells to cut and pattern of cuts
00462         cellCuts
00463         (
00464             const polyMesh& mesh,
00465             const labelList& cutCells,
00466             const labelList& meshVerts,
00467             const labelList& meshEdges,
00468             const scalarField& meshEdgeWeights
00469         );
00470 
00471         //- Construct from pattern of cuts. Detect cells to cut.
00472         cellCuts
00473         (
00474             const polyMesh& mesh,
00475             const labelList& meshVerts,
00476             const labelList& meshEdges,
00477             const scalarField& meshEdgeWeights
00478         );
00479 
00480         //- Construct from complete cellLoops through specified cells.
00481         //  Checks for consistency.
00482         //  Constructs cut-cut addressing and cellAnchorPoints.
00483         cellCuts
00484         (
00485             const polyMesh& mesh,
00486             const labelList& cellLabels,
00487             const labelListList& cellLoops,
00488             const List<scalarField>& cellEdgeWeights
00489         );
00490 
00491         //- Construct from list of cells to cut and direction to cut in
00492         //  (always through cell centre) and celllooper.
00493         cellCuts
00494         (
00495             const polyMesh& mesh,
00496             const cellLooper& cellCutter,
00497             const List<refineCell>& refCells
00498         );
00499 
00500         //- Construct from list of cells to cut and plane to cut with and
00501         //  celllooper. (constructor above always cuts through cell centre)
00502         cellCuts
00503         (
00504             const polyMesh& mesh,
00505             const cellLooper& cellCutter,
00506             const labelList& cellLabels,
00507             const List<plane>& cutPlanes
00508         );
00509 
00510         //- Construct from components
00511         cellCuts
00512         (
00513             const polyMesh& mesh,
00514             const boolList& pointIsCut,
00515             const boolList& edgeIsCut,
00516             const scalarField& edgeWeight,
00517             const Map<edge>& faceSplitCut,
00518             const labelListList& cellLoops,
00519             const label nLoops,
00520             const labelListList& cellAnchorPoints
00521         );
00522 
00523 
00524     // Destructor
00525 
00526         ~cellCuts();
00527 
00528         //- Clear out demand driven storage
00529         void clearOut();
00530 
00531 
00532     // Member Functions
00533 
00534         // Access
00535 
00536             //- Is mesh point cut
00537             const boolList& pointIsCut() const
00538             {
00539                 return pointIsCut_;
00540             }
00541 
00542             //- Is edge cut
00543             const boolList& edgeIsCut() const
00544             {
00545                 return edgeIsCut_;
00546             }
00547 
00548             //- If edge is cut gives weight (ratio between start() and end())
00549             const scalarField& edgeWeight() const
00550             {
00551                 return edgeWeight_;
00552             }
00553 
00554             //- Cuts per existing face (includes those along edge of face)
00555             //  Cuts in no particular order
00556             const labelListList& faceCuts() const
00557             {
00558                 if (!faceCutsPtr_)
00559                 {
00560                     calcFaceCuts();
00561                 }
00562                 return *faceCutsPtr_;
00563             }
00564 
00565             //- Gives for split face the two cuts that split the face into two.
00566             const Map<edge>& faceSplitCut() const
00567             {
00568                 return faceSplitCut_;
00569             }
00570 
00571             //- For each cut cell the cut along the circumference.
00572             const labelListList& cellLoops() const
00573             {
00574                 return cellLoops_;
00575             }
00576 
00577             //- Number of valid cell loops
00578             label nLoops() const
00579             {
00580                 return nLoops_;
00581             }
00582 
00583             //- For each cut cell the points on the 'anchor' side of the cell.
00584             const labelListList& cellAnchorPoints() const
00585             {
00586                 return cellAnchorPoints_;
00587             }
00588 
00589 
00590         // Other
00591 
00592             //- Returns coordinates of points on loop for given cell.
00593             //  Uses cellLoops_ and edgeWeight_
00594             pointField loopPoints(const label cellI) const;
00595   
00596             //- Invert anchor point selection.
00597             labelList nonAnchorPoints
00598             (
00599                 const labelList& cellPoints,
00600                 const labelList& anchorPoints,
00601                 const labelList& loop
00602             ) const;
00603 
00604             //- Flip loop for cellI. Updates anchor points as well.
00605             void flip(const label cellI);
00606 
00607             //- Flip loop for cellI. Does not update anchors. Use with care
00608             //  (only if you're sure loop orientation is wrong)
00609             void flipLoopOnly(const label cellI);
00610 
00611 
00612         // Write
00613 
00614             //- debugging:Write list of cuts to stream in OBJ format
00615             void writeOBJ
00616             (
00617                 Ostream& os,
00618                 const pointField& loopPoints,
00619                 label& vertI
00620             ) const;
00621 
00622             //- debugging:Write all of cuts to stream in OBJ format
00623             void writeOBJ(Ostream& os) const;
00624 
00625             //- debugging:Write edges of cell and loop
00626             void writeCellOBJ(const fileName& dir, const label cellI) const;
00627 
00628 };
00629 
00630 
00631 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00632 
00633 } // End namespace Foam
00634 
00635 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00636 
00637 #endif
00638 
00639 // ************************************************************************* //
For further information go to www.openfoam.org