applications/utilities/postProcessing/dataConversion/foamToFieldview9/fieldviewTopology.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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 \*---------------------------------------------------------------------------*/
  25
  26 #include "fieldviewTopology.H"
  27 #include "polyMesh.H"
  28 #include "cellShape.H"
  29 #include "cellModeller.H"
  30 #include "wallPolyPatch.H"
  31 #include "symmetryPolyPatch.H"
  32
  33
  34 #include "fv_reader_tags.h"
  35
  36 extern "C"
  37 {
  38     unsigned int fv_encode_elem_header(int elem_type, int wall_info[]);
  39 }
  40
  41
  42 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  43
  44 Foam::labelList Foam::fieldviewTopology::calcFaceAddressing
  45 (
  46     const faceList& allFaces,   // faces given faceLabels
  47     const cellShape& shape,
  48     const labelList& faces,     // faceLabels for given cell
  49     const label cellI
  50 )
  51 {
  52     // return value.
  53     labelList shapeToMesh(shape.nFaces(), -1);
  54
  55     const faceList modelFaces(shape.faces());
  56
  57     // Loop over all faces of cellShape
  58     forAll(modelFaces, cellFaceI)
  59     {
  60         // face (vertex list)
  61         const face& modelFace = modelFaces[cellFaceI];
  62
  63         // Loop over all face labels
  64         forAll(faces, faceI)
  65         {
  66             const face& vertLabels = allFaces[faces[faceI]];
  67
  68             if (vertLabels == modelFace)
  69             {
  70                 shapeToMesh[cellFaceI] = faces[faceI];
  71                 break;
  72             }
  73         }
  74
  75         if (shapeToMesh[cellFaceI] == -1)
  76         {
  77             FatalErrorIn("foamToFieldview : calcFaceAddressing")
  78                 << "calcFaceAddressing : can't match face to shape.\n"
  79                 << "    shape face:" << modelFace << endl
  80                 << "    face labels:" << faces << endl
  81                 << "    cellI:" << cellI << endl;
  82
  83             FatalError << "Faces consist of vertices:" << endl;
  84             forAll(faces, faceI)
  85             {
  86                 FatalError
  87                     << "    face:" << faces[faceI]
  88                     << allFaces[faces[faceI]] << endl;
  89             }
  90             FatalError << exit(FatalError);
  91         }
  92     }
  93     return shapeToMesh;
  94 }
  95
  96
  97 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  98
  99 // Construct from components
 100 Foam::fieldviewTopology::fieldviewTopology
 101 (
 102     const polyMesh& mesh,
 103     const bool setWallInfo
 104 )
 105 :
 106     hexLabels_((1+8)*mesh.nCells()),
 107     prismLabels_((1+6)*mesh.nCells()),
 108     pyrLabels_((1+5)*mesh.nCells()),
 109     tetLabels_((1+4)*mesh.nCells()),
 110     nPoly_(0),
 111     quadFaceLabels_(mesh.boundaryMesh().size()),
 112     nPolyFaces_(mesh.boundaryMesh().size())
 113 {
 114     // Mark all faces that are to be seen as wall for particle
 115     // tracking and all cells that use one or more of these walls
 116
 117     labelList wallFace(mesh.nFaces(), NOT_A_WALL);
 118     boolList wallCell(mesh.nCells(), false);
 119
 120     if (setWallInfo)
 121     {
 122         forAll(mesh.boundaryMesh(), patchI)
 123         {
 124             const polyPatch& currPatch = mesh.boundaryMesh()[patchI];
 125             if
 126             (
 127                 isA<wallPolyPatch>(currPatch)
 128              || isA<symmetryPolyPatch>(currPatch)
 129             )
 130             {
 131                 forAll(currPatch, patchFaceI)
 132                 {
 133                     label meshFaceI = currPatch.start() + patchFaceI;
 134
 135                     wallFace[meshFaceI] = A_WALL;
 136                     wallCell[mesh.faceOwner()[meshFaceI]] = true;
 137                 }
 138             }
 139         }
 140     }
 141
 142
 143
 144     const cellModel& tet = *(cellModeller::lookup("tet"));
 145     const cellModel& pyr = *(cellModeller::lookup("pyr"));
 146     const cellModel& prism = *(cellModeller::lookup("prism"));
 147     const cellModel& wedge = *(cellModeller::lookup("wedge"));
 148     const cellModel& tetWedge = *(cellModeller::lookup("tetWedge"));
 149     const cellModel& hex = *(cellModeller::lookup("hex"));
 150
 151     // Pre calculate headers for cells not on walls
 152     labelList notWallFlags(6, NOT_A_WALL);
 153     label tetNotWall = fv_encode_elem_header
 154     (
 155         FV_TET_ELEM_ID, notWallFlags.begin()
 156     );
 157     label pyrNotWall = fv_encode_elem_header
 158     (
 159         FV_PYRA_ELEM_ID, notWallFlags.begin()
 160     );
 161     label prismNotWall = fv_encode_elem_header
 162     (
 163         FV_PRISM_ELEM_ID, notWallFlags.begin()
 164     );
 165     label hexNotWall = fv_encode_elem_header
 166     (
 167         FV_HEX_ELEM_ID, notWallFlags.begin()
 168     );
 169
 170     // Some aliases
 171     const cellList& cellFaces = mesh.cells();
 172     const cellShapeList& cellShapes = mesh.cellShapes();
 173
 174
 175     label hexi = 0;
 176     label prismi = 0;
 177     label pyri = 0;
 178     label teti = 0;
 179
 180     const faceList& allFaces = mesh.faces();
 181
 182     labelList wallFlags(6);
 183     forAll(cellShapes, celli)
 184     {
 185         const cellShape& cellShape = cellShapes[celli];
 186         const cellModel& cellModel = cellShape.model();
 187
 188         if (cellModel == tet)
 189         {
 190             if (!wallCell[celli])
 191             {
 192                 tetLabels_[teti++] = tetNotWall;
 193             }
 194             else
 195             {
 196                 labelList modelToMesh = calcFaceAddressing
 197                 (
 198                     allFaces, cellShape, cellFaces[celli], celli
 199                 );
 200
 201                 wallFlags[0] = wallFace[modelToMesh[0]];
 202                 wallFlags[1] = wallFace[modelToMesh[1]];
 203                 wallFlags[2] = wallFace[modelToMesh[2]];
 204                 wallFlags[3] = wallFace[modelToMesh[3]];
 205
 206                 tetLabels_[teti++] = fv_encode_elem_header
 207                 (
 208                     FV_TET_ELEM_ID, wallFlags.begin()
 209                 );
 210             }
 211
 212             tetLabels_[teti++] = cellShape[0] + 1;
 213             tetLabels_[teti++] = cellShape[1] + 1;
 214             tetLabels_[teti++] = cellShape[2] + 1;
 215             tetLabels_[teti++] = cellShape[3] + 1;
 216         }
 217         else if (cellModel == pyr)
 218         {
 219             if (!wallCell[celli])
 220             {
 221                 pyrLabels_[pyri++] = pyrNotWall;
 222             }
 223             else
 224             {
 225                 labelList modelToMesh = calcFaceAddressing
 226                 (
 227                     allFaces, cellShape, cellFaces[celli], celli
 228                 );
 229
 230                 wallFlags[0] = wallFace[modelToMesh[0]];
 231                 wallFlags[1] = wallFace[modelToMesh[3]];
 232                 wallFlags[2] = wallFace[modelToMesh[2]];
 233                 wallFlags[3] = wallFace[modelToMesh[1]];
 234                 wallFlags[4] = wallFace[modelToMesh[4]];
 235
 236                 pyrLabels_[pyri++] = fv_encode_elem_header
 237                 (
 238                     FV_PYRA_ELEM_ID, wallFlags.begin()
 239                 );
 240             }
 241
 242             pyrLabels_[pyri++] = cellShape[0] + 1;
 243             pyrLabels_[pyri++] = cellShape[1] + 1;
 244             pyrLabels_[pyri++] = cellShape[2] + 1;
 245             pyrLabels_[pyri++] = cellShape[3] + 1;
 246             pyrLabels_[pyri++] = cellShape[4] + 1;
 247         }
 248         else if (cellModel == prism)
 249         {
 250             if (!wallCell[celli])
 251             {
 252                 prismLabels_[prismi++] = prismNotWall;
 253             }
 254             else
 255             {
 256                 labelList modelToMesh = calcFaceAddressing
 257                 (
 258                     allFaces, cellShape, cellFaces[celli], celli
 259                 );
 260
 261                 wallFlags[0] = wallFace[modelToMesh[4]];
 262                 wallFlags[1] = wallFace[modelToMesh[2]];
 263                 wallFlags[2] = wallFace[modelToMesh[3]];
 264                 wallFlags[3] = wallFace[modelToMesh[0]];
 265                 wallFlags[4] = wallFace[modelToMesh[1]];
 266
 267                 prismLabels_[prismi++] = fv_encode_elem_header
 268                 (
 269                     FV_PRISM_ELEM_ID, wallFlags.begin()
 270                 );
 271             }
 272
 273             prismLabels_[prismi++] = cellShape[0] + 1;
 274             prismLabels_[prismi++] = cellShape[3] + 1;
 275             prismLabels_[prismi++] = cellShape[4] + 1;
 276             prismLabels_[prismi++] = cellShape[1] + 1;
 277             prismLabels_[prismi++] = cellShape[5] + 1;
 278             prismLabels_[prismi++] = cellShape[2] + 1;
 279         }
 280         else if (cellModel == tetWedge)
 281         {
 282             // Treat as prism with collapsed edge
 283             if (!wallCell[celli])
 284             {
 285                 prismLabels_[prismi++] = prismNotWall;
 286             }
 287             else
 288             {
 289                 labelList modelToMesh = calcFaceAddressing
 290                 (
 291                     allFaces, cellShape, cellFaces[celli], celli
 292                 );
 293
 294                 wallFlags[0] = wallFace[modelToMesh[1]];
 295                 wallFlags[1] = wallFace[modelToMesh[2]];
 296                 wallFlags[2] = wallFace[modelToMesh[3]];
 297                 wallFlags[3] = wallFace[modelToMesh[0]];
 298                 wallFlags[4] = wallFace[modelToMesh[3]];
 299
 300                 prismLabels_[prismi++] = fv_encode_elem_header
 301                 (
 302                     FV_PRISM_ELEM_ID, wallFlags.begin()
 303                 );
 304             }
 305
 306             prismLabels_[prismi++] = cellShape[0] + 1;
 307             prismLabels_[prismi++] = cellShape[3] + 1;
 308             prismLabels_[prismi++] = cellShape[4] + 1;
 309             prismLabels_[prismi++] = cellShape[1] + 1;
 310             prismLabels_[prismi++] = cellShape[4] + 1;
 311             prismLabels_[prismi++] = cellShape[2] + 1;
 312         }
 313         else if (cellModel == wedge)
 314         {
 315             if (!wallCell[celli])
 316             {
 317                 hexLabels_[hexi++] = hexNotWall;
 318             }
 319             else
 320             {
 321                 labelList modelToMesh = calcFaceAddressing
 322                 (
 323                     allFaces, cellShape, cellFaces[celli], celli
 324                 );
 325
 326                 wallFlags[0] = wallFace[modelToMesh[2]];
 327                 wallFlags[1] = wallFace[modelToMesh[3]];
 328                 wallFlags[2] = wallFace[modelToMesh[0]];
 329                 wallFlags[3] = wallFace[modelToMesh[1]];
 330                 wallFlags[4] = wallFace[modelToMesh[4]];
 331                 wallFlags[5] = wallFace[modelToMesh[5]];
 332
 333                 hexLabels_[hexi++] = fv_encode_elem_header
 334                 (
 335                     FV_HEX_ELEM_ID, wallFlags.begin()
 336                 );
 337             }
 338             hexLabels_[hexi++] = cellShape[0] + 1;
 339             hexLabels_[hexi++] = cellShape[1] + 1;
 340             hexLabels_[hexi++] = cellShape[0] + 1;
 341             hexLabels_[hexi++] = cellShape[2] + 1;
 342             hexLabels_[hexi++] = cellShape[3] + 1;
 343             hexLabels_[hexi++] = cellShape[4] + 1;
 344             hexLabels_[hexi++] = cellShape[6] + 1;
 345             hexLabels_[hexi++] = cellShape[5] + 1;
 346         }
 347         else if (cellModel == hex)
 348         {
 349             if (!wallCell[celli])
 350             {
 351                 hexLabels_[hexi++] = hexNotWall;
 352             }
 353             else
 354             {
 355                 labelList modelToMesh = calcFaceAddressing
 356                 (
 357                     allFaces, cellShape, cellFaces[celli], celli
 358                 );
 359
 360                 wallFlags[0] = wallFace[modelToMesh[0]];
 361                 wallFlags[1] = wallFace[modelToMesh[1]];
 362                 wallFlags[2] = wallFace[modelToMesh[4]];
 363                 wallFlags[3] = wallFace[modelToMesh[5]];
 364                 wallFlags[4] = wallFace[modelToMesh[2]];
 365                 wallFlags[5] = wallFace[modelToMesh[3]];
 366
 367                 hexLabels_[hexi++] = fv_encode_elem_header
 368                 (
 369                     FV_HEX_ELEM_ID, wallFlags.begin()
 370                 );
 371             }
 372             hexLabels_[hexi++] = cellShape[0] + 1;
 373             hexLabels_[hexi++] = cellShape[1] + 1;
 374             hexLabels_[hexi++] = cellShape[3] + 1;
 375             hexLabels_[hexi++] = cellShape[2] + 1;
 376             hexLabels_[hexi++] = cellShape[4] + 1;
 377             hexLabels_[hexi++] = cellShape[5] + 1;
 378             hexLabels_[hexi++] = cellShape[7] + 1;
 379             hexLabels_[hexi++] = cellShape[6] + 1;
 380         }
 381         else
 382         {
 383             nPoly_++;
 384         }
 385     }
 386
 387     hexLabels_.setSize(hexi);
 388     prismLabels_.setSize(prismi);
 389     pyrLabels_.setSize(pyri);
 390     tetLabels_.setSize(teti);
 391
 392
 393     //
 394     // Patches
 395     //
 396     forAll(mesh.boundaryMesh(), patchI)
 397     {
 398         const polyPatch& patchFaces = mesh.boundaryMesh()[patchI];
 399
 400         labelList& faceLabels = quadFaceLabels_[patchI];
 401
 402         // Faces, each 4 labels. Size big enough
 403         faceLabels.setSize(patchFaces.size()*4);
 404
 405         label labelI = 0;
 406
 407         forAll(patchFaces, faceI)
 408         {
 409             const face& patchFace = patchFaces[faceI];
 410
 411             if (patchFace.size() == 3)
 412             {
 413                 faceLabels[labelI++] = patchFace[0] + 1;
 414                 faceLabels[labelI++] = patchFace[1] + 1;
 415                 faceLabels[labelI++] = patchFace[2] + 1;
 416                 faceLabels[labelI++] = 0;   // Fieldview:triangle definition
 417             }
 418             else if (patchFace.size() == 4)
 419             {
 420                 faceLabels[labelI++] = patchFace[0] + 1;
 421                 faceLabels[labelI++] = patchFace[1] + 1;
 422                 faceLabels[labelI++] = patchFace[2] + 1;
 423                 faceLabels[labelI++] = patchFace[3] + 1;
 424             }
 425         }
 426
 427         faceLabels.setSize(labelI);
 428
 429         label nFaces = labelI/4;
 430
 431         nPolyFaces_[patchI] = patchFaces.size() - nFaces;
 432     }
 433 }
 434
 435
 436 // ************************************************************************* //