src/sampling/sampledSurface/isoSurface/isoSurfaceCell.H

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