applications/utilities/postProcessing/dataConversion/foamToVTK/vtkTopo.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2010 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
  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 Description
  25
  26     Note: bug in vtk displaying wedges? Seems to display ok if we decompose
  27     them. Should be thoroughly tested!
  28     (they appear rarely in polyhedral meshes, do appear in some cut meshes)
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 #include "vtkTopo.H"
  33 #include <OpenFOAM/polyMesh.H>
  34 #include <OpenFOAM/cellShape.H>
  35 #include <OpenFOAM/cellModeller.H>
  36
  37 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  38
  39 // Construct from components
  40 Foam::vtkTopo::vtkTopo(const polyMesh& mesh)
  41 :
  42     mesh_(mesh),
  43     vertLabels_(),
  44     cellTypes_(),
  45     addPointCellLabels_(),
  46     superCells_()
  47 {
  48     const cellModel& tet = *(cellModeller::lookup("tet"));
  49     const cellModel& pyr = *(cellModeller::lookup("pyr"));
  50     const cellModel& prism = *(cellModeller::lookup("prism"));
  51     const cellModel& tetWedge = *(cellModeller::lookup("tetWedge"));
  52     const cellModel& hex = *(cellModeller::lookup("hex"));
  53
  54
  55     const cellShapeList& cellShapes = mesh_.cellShapes();
  56
  57
  58     // Number of additional points needed by the decomposition of polyhedra
  59     label nAddPoints = 0;
  60
  61     // Number of additional cells generated by the decomposition of polyhedra
  62     label nAddCells = 0;
  63
  64     // Scan for cells which need to be decomposed and count additional points
  65     // and cells
  66
  67     forAll(cellShapes, cellI)
  68     {
  69         const cellModel& model = cellShapes[cellI].model();
  70
  71         if
  72         (
  73             model != hex
  74 //         && model != wedge            // See above.
  75          && model != prism
  76          && model != pyr
  77          && model != tet
  78          && model != tetWedge
  79         )
  80         {
  81             const cell& cFaces = mesh_.cells()[cellI];
  82
  83             forAll(cFaces, cFaceI)
  84             {
  85                 const face& f = mesh_.faces()[cFaces[cFaceI]];
  86
  87                 label nQuads = 0;
  88                 label nTris = 0;
  89                 f.nTrianglesQuads(mesh_.points(), nTris, nQuads);
  90
  91                 nAddCells += nQuads + nTris;
  92             }
  93
  94             nAddCells--;
  95             nAddPoints++;
  96         }
  97     }
  98
  99     // Set size of additional point addressing array
 100     // (from added point to original cell)
 101     addPointCellLabels_.setSize(nAddPoints);
 102
 103     // Set size of additional cells mapping array
 104     // (from added cell to original cell)
 105     superCells_.setSize(nAddCells);
 106
 107     // List of vertex labels in VTK ordering
 108     vertLabels_.setSize(cellShapes.size() + nAddCells);
 109
 110     // Label of vtk type
 111     cellTypes_.setSize(cellShapes.size() + nAddCells);
 112
 113     // Set counters for additional points and additional cells
 114     label api = 0, aci = 0;
 115
 116     forAll(cellShapes, cellI)
 117     {
 118         const cellShape& cellShape = cellShapes[cellI];
 119         const cellModel& cellModel = cellShape.model();
 120
 121         labelList& vtkVerts = vertLabels_[cellI];
 122
 123         if (cellModel == tet)
 124         {
 125             vtkVerts = cellShape;
 126
 127             cellTypes_[cellI] = VTK_TETRA;
 128         }
 129         else if (cellModel == pyr)
 130         {
 131             vtkVerts = cellShape;
 132
 133             cellTypes_[cellI] = VTK_PYRAMID;
 134         }
 135         else if (cellModel == prism)
 136         {
 137             vtkVerts.setSize(6);
 138             vtkVerts[0] = cellShape[0];
 139             vtkVerts[1] = cellShape[2];
 140             vtkVerts[2] = cellShape[1];
 141             vtkVerts[3] = cellShape[3];
 142             vtkVerts[4] = cellShape[5];
 143             vtkVerts[5] = cellShape[4];
 144
 145             // VTK calls this a wedge.
 146             cellTypes_[cellI] = VTK_WEDGE;
 147         }
 148         else if (cellModel == tetWedge)
 149         {
 150             // Treat as squeezed prism
 151             vtkVerts.setSize(6);
 152             vtkVerts[0] = cellShape[0];
 153             vtkVerts[1] = cellShape[2];
 154             vtkVerts[2] = cellShape[1];
 155             vtkVerts[3] = cellShape[3];
 156             vtkVerts[4] = cellShape[4];
 157             vtkVerts[5] = cellShape[4];
 158
 159             cellTypes_[cellI] = VTK_WEDGE;
 160         }
 161 //        else if (cellModel == wedge)
 162 //        {
 163 //            // Treat as squeezed hex
 164 //            vtkVerts.setSize(8);
 165 //            vtkVerts[0] = cellShape[0];
 166 //            vtkVerts[1] = cellShape[1];
 167 //            vtkVerts[2] = cellShape[2];
 168 //            vtkVerts[3] = cellShape[0];
 169 //            vtkVerts[4] = cellShape[3];
 170 //            vtkVerts[5] = cellShape[4];
 171 //            vtkVerts[6] = cellShape[5];
 172 //            vtkVerts[7] = cellShape[6];
 173 //
 174 //            cellTypes_[cellI] = VTK_HEXAHEDRON;
 175 //        }
 176         else if (cellModel == hex)
 177         {
 178             vtkVerts.setSize(8);
 179             vtkVerts[0] = cellShape[0];
 180             vtkVerts[1] = cellShape[1];
 181             vtkVerts[2] = cellShape[2];
 182             vtkVerts[3] = cellShape[3];
 183             vtkVerts[4] = cellShape[4];
 184             vtkVerts[5] = cellShape[5];
 185             vtkVerts[6] = cellShape[6];
 186             vtkVerts[7] = cellShape[7];
 187
 188             cellTypes_[cellI] = VTK_HEXAHEDRON;
 189         }
 190         else
 191         {
 192             // Polyhedral cell. Decompose into tets + prisms.
 193             // (see dxFoamExec/createDxConnections.C)
 194
 195             // Mapping from additional point to cell
 196             addPointCellLabels_[api] = cellI;
 197
 198             // Whether to insert cell in place of original or not.
 199             bool substituteCell = true;
 200
 201             const labelList& cFaces = mesh_.cells()[cellI];
 202
 203             forAll(cFaces, cFaceI)
 204             {
 205                 const face& f = mesh_.faces()[cFaces[cFaceI]];
 206
 207                 // Number of triangles and quads in decomposition
 208                 label nTris = 0;
 209                 label nQuads = 0;
 210                 f.nTrianglesQuads(mesh_.points(), nTris, nQuads);
 211
 212                 // Do actual decomposition into triFcs and quadFcs.
 213                 faceList triFcs(nTris);
 214                 faceList quadFcs(nQuads);
 215                 label trii = 0;
 216                 label quadi = 0;
 217                 f.trianglesQuads(mesh_.points(), trii, quadi, triFcs, quadFcs);
 218
 219                 forAll(quadFcs, quadi)
 220                 {
 221                     label thisCellI = -1;
 222
 223                     if (substituteCell)
 224                     {
 225                         thisCellI = cellI;
 226
 227                         substituteCell = false;
 228                     }
 229                     else
 230                     {
 231                         thisCellI = mesh_.nCells() + aci;
 232
 233                         superCells_[aci] = cellI;
 234
 235                         aci++;
 236                     }
 237
 238                     labelList& addVtkVerts = vertLabels_[thisCellI];
 239
 240                     addVtkVerts.setSize(5);
 241
 242                     const face& quad = quadFcs[quadi];
 243
 244                     addVtkVerts[0] = quad[0];
 245                     addVtkVerts[1] = quad[1];
 246                     addVtkVerts[2] = quad[2];
 247                     addVtkVerts[3] = quad[3];
 248                     addVtkVerts[4] = mesh_.nPoints() + api;
 249
 250                     cellTypes_[thisCellI] = VTK_PYRAMID;
 251                 }
 252
 253                 forAll(triFcs, trii)
 254                 {
 255                     label thisCellI = -1;
 256
 257                     if (substituteCell)
 258                     {
 259                         thisCellI = cellI;
 260
 261                         substituteCell = false;
 262                     }
 263                     else
 264                     {
 265                         thisCellI = mesh_.nCells() + aci;
 266
 267                         superCells_[aci] = cellI;
 268
 269                         aci++;
 270                     }
 271
 272
 273                     labelList& addVtkVerts = vertLabels_[thisCellI];
 274
 275                     const face& tri = triFcs[trii];
 276
 277                     addVtkVerts.setSize(4);
 278                     addVtkVerts[0] = tri[0];
 279                     addVtkVerts[1] = tri[1];
 280                     addVtkVerts[2] = tri[2];
 281                     addVtkVerts[3] = mesh_.nPoints() + api;
 282
 283                     cellTypes_[thisCellI] = VTK_TETRA;
 284                 }
 285             }
 286
 287             api++;
 288         }
 289     }
 290
 291     Pout<< "    Original cells:" << mesh_.nCells()
 292         << " points:" << mesh_.nPoints()
 293         << "   Additional cells:" << superCells_.size()
 294         << "  additional points:" << addPointCellLabels_.size()
 295         << nl << endl;
 296 }
 297
 298 // ************************ vim: set sw=4 sts=4 et: ************************ //