src/sampling/sampledSurface/isoSurface/isoSurfaceCell.H

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend 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 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  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.H"
  48 #include "labelPair.H"
  49 #include "pointIndexHit.H"
  50 #include "PackedBoolList.H"
  51
  52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  53
  54 namespace Foam
  55 {
  56
  57 class polyMesh;
  58
  59 /*---------------------------------------------------------------------------*\
  60                        Class isoSurfaceCell Declaration
  61 \*---------------------------------------------------------------------------*/
  62
  63 class isoSurfaceCell
  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         (
 143             pointField& localPoints,
 144             DynamicList<labelledTri, 64>& localTris
 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,
 219             DynamicList<Type>& points
 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             dynamicLabelList& 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,
 254             List<FixedList<label, 3> >& faceEdges,
 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,
 264             const List<FixedList<label, 3> >& faceEdges,
 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&,
 275             const List<FixedList<label, 3> >& faceEdges,
 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         (
 291             const List<FixedList<label, 3> >& faceEdges,
 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
 318         isoSurfaceCell
 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
 338         const labelList triPointMergeMap() const
 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>
 347         tmp<Field<Type> > interpolate
 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 // ************************************************************************* //