applications/utilities/mesh/conversion/ansysToFoam/ansysToFoam.L

   1 /*--------------------------------*- C++ -*----------------------------------*\
   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 Application
  25     ansysToFoam
  26
  27 Description
  28     Converts an ANSYS input mesh file, exported from I-DEAS, to FOAM format.
  29
  30 \*---------------------------------------------------------------------------*/
  31
  32 %{
  33
  34 #undef yyFlexLexer
  35
  36  /* ------------------------------------------------------------------------ *\
  37    ------ local definitions
  38  \* ------------------------------------------------------------------------ */
  39
  40 #include <sstream>
  41 // For EOF only
  42 #include <cstdio>
  43
  44 #include "scalar.H"
  45 #include "IStringStream.H"
  46
  47 using namespace Foam;
  48
  49 #include "argList.H"
  50 #include "objectRegistry.H"
  51 #include "foamTime.H"
  52 #include "polyMesh.H"
  53 #include "emptyPolyPatch.H"
  54 #include "preservePatchTypes.H"
  55 #include "cellShape.H"
  56 #include "cellModeller.H"
  57 #include "SLList.H"
  58 #include "SLPtrList.H"
  59
  60 SLList<point> slPoints;
  61 SLList<label> slPointMap;
  62 label maxNodei = 0;
  63
  64 SLPtrList<labelList> slCellLabels;
  65 SLList<label> slCellMap;
  66 SLList<label> slCellType;
  67 label maxCelli = 0;
  68
  69 PtrList<SLList<label> > slPatchCells;
  70 PtrList<SLList<label> > slPatchCellFaces;
  71
  72 // Cell types
  73 Map<word> cellTypes;
  74 label currentTypei = -1;
  75
  76
  77 // Dummy yywrap to keep yylex happy at compile time.
  78 // It is called by yylex but is not used as the mechanism to change file.
  79 // See <<EOF>>
  80 #if YY_FLEX_MINOR_VERSION < 6 && YY_FLEX_SUBMINOR_VERSION < 34
  81 extern "C" int yywrap()
  82 #else
  83 int yyFlexLexer::yywrap()
  84 #endif
  85 {
  86     return 1;
  87 }
  88
  89 %}
  90
  91 one_space             [ \t\f\r]
  92 space                 {one_space}*
  93 some_space            {one_space}+
  94 cspace                ","{space}
  95
  96 alpha                 [_A-Za-z]
  97 digit                 [0-9]
  98
  99 identifier            {alpha}({alpha}|{digit})*
 100 integer               {digit}+
 101 label                 [1-9]{digit}*
 102
 103 exponent_part         [eE][-+]?{digit}+
 104 fractional_constant   [-+]?(({digit}*"."{digit}+)|({digit}+"."?))
 105
 106 floatNum              (({fractional_constant}{exponent_part}?)|({digit}+{exponent_part}))
 107
 108 x                     {floatNum}
 109 y                     {floatNum}
 110 z                     {floatNum}
 111 value                 {floatNum}
 112
 113 node                  ^{space}"N"{cspace}
 114 element               ^{space}"EN"{cspace}
 115 bface                 ^{space}"SFE"{cspace}
 116 elementTypeName       ^{space}"ET"{cspace}
 117 elementType           ^{space}"TYPE"{cspace}
 118
 119
 120 %%
 121
 122 %{
 123     labelList labels(8);
 124 %}
 125
 126
 127  /* ------------------------------------------------------------------------ *\
 128                             ------ Start Lexing ------
 129  \* ------------------------------------------------------------------------ */
 130
 131 {node}{label}{cspace}{x}{cspace}{y}{cspace}{z}{space}\n {
 132         IStringStream nodeStream(YYText());
 133         char tag, c;
 134         label nodei;
 135         point node;
 136         nodeStream
 137             >> tag
 138             >> c >> nodei
 139             >> c >> node.x()
 140             >> c >> node.y()
 141             >> c >> node.z();
 142
 143         if (nodei > maxNodei) maxNodei = nodei;
 144
 145         slPointMap.append(nodei);
 146         slPoints.append(node);
 147     }
 148
 149
 150 {element}{label}{cspace}{label}{cspace}{label}{cspace}{label}{cspace}{label}{cspace}{label}{cspace}{label}{cspace}{label}{cspace}{label}{space}\n {
 151         IStringStream elementStream(YYText());
 152         char tag, c;
 153         label celli;
 154         elementStream
 155             >> tag >> tag
 156             >> c >> celli
 157             >> c >> labels[0]
 158             >> c >> labels[1]
 159             >> c >> labels[2]
 160             >> c >> labels[3]
 161             >> c >> labels[4]
 162             >> c >> labels[5]
 163             >> c >> labels[6]
 164             >> c >> labels[7];
 165
 166         if (celli > maxCelli) maxCelli = celli;
 167
 168         slCellMap.append(celli);
 169         slCellLabels.append(new labelList(labels));
 170         slCellType.append(currentTypei);
 171     }
 172
 173
 174 {bface}{label}{cspace}{label}{cspace}{identifier}{cspace}{integer}{cspace}{value}{space}\n {
 175         IStringStream bfaceStream(YYText());
 176         char tag, c;
 177         label elementi;
 178         label facei;
 179         scalar indexValue, unknown;
 180         bfaceStream
 181             >> tag >> tag >> tag
 182             >> c >> elementi
 183             >> c >> facei
 184             >> c >> tag >> tag >> tag >> tag
 185             >> c >> unknown
 186             >> c >> indexValue;
 187
 188         label patchi = label(indexValue);
 189
 190         if (patchi > slPatchCells.size())
 191         {
 192             slPatchCells.setSize(patchi);
 193
 194             forAll(slPatchCells, i)
 195             {
 196                 if (!slPatchCells(i))
 197                 {
 198                     slPatchCells.set(i, new SLList<label>);
 199                 }
 200             }
 201         }
 202
 203         if (patchi > slPatchCellFaces.size())
 204         {
 205             slPatchCellFaces.setSize(patchi);
 206
 207             forAll(slPatchCells, i)
 208             {
 209                 if (!slPatchCellFaces(i))
 210                 {
 211                     slPatchCellFaces.set(i, new SLList<label>);
 212                 }
 213             }
 214         }
 215
 216         slPatchCells[patchi-1].append(elementi);
 217         slPatchCellFaces[patchi-1].append(facei);
 218     }
 219
 220
 221 {elementTypeName}{label}{cspace}{identifier}{space}\n {
 222
 223         IStringStream elementStream(YYText());
 224         char tag,c;
 225         label cellTypei;
 226         word cellTypeName;
 227         elementStream
 228             >> tag >> tag           // skip 'ET'
 229             >> c >> cellTypei
 230             >> c >> cellTypeName;
 231
 232         Info<< "Read typeName " << cellTypeName
 233             << " for type " << cellTypei << endl;
 234
 235         cellTypes.insert(cellTypei, cellTypeName);
 236     }
 237
 238
 239 {elementType}{label}{space}\n {
 240         IStringStream elementStream(YYText());
 241         char tag,c;
 242         label cellTypei;
 243         elementStream
 244             >> tag >> tag >> tag >> tag     // skip 'TYPE'
 245             >> c >> cellTypei;
 246
 247         currentTypei = cellTypei;
 248     }
 249
 250
 251
 252  /* ------------------------------------------------------------------------ *\
 253     ------ Ignore remaining space and \n s.  Any other characters are errors.
 254  \* ------------------------------------------------------------------------ */
 255
 256 .|\n {}
 257
 258
 259  /* ------------------------------------------------------------------------ *\
 260     ------ On EOF return to previous file, if none exists terminate.
 261  \* ------------------------------------------------------------------------ */
 262
 263 <<EOF>> {
 264             yyterminate();
 265     }
 266 %%
 267
 268
 269 #include "fileName.H"
 270 #include <fstream>
 271 using std::ifstream;
 272
 273 label findFace(const polyMesh& mesh, const face& f)
 274 {
 275     if (f.empty())
 276     {
 277         FatalErrorIn("findFace(const polyMesh&, const face&)")
 278             << "Empty face found: " << f
 279             << exit(FatalError);
 280     }
 281
 282     const faceList& faces = mesh.faces();
 283
 284     const labelList& pFaces = mesh.pointFaces()[f[0]];
 285
 286     forAll (pFaces, i)
 287     {
 288         label faceI = pFaces[i];
 289
 290         if (faces[faceI] == f)
 291         {
 292             return faceI;
 293         }
 294     }
 295
 296     FatalErrorIn("findFace(const polyMesh&, const face&)")
 297         << "Cannot find a face matching " << f
 298         << exit(FatalError);
 299
 300     return -1;
 301 }
 302
 303
 304 int main(int argc, char *argv[])
 305 {
 306     argList::noParallel();
 307     argList::validArgs.append("ANSYS input file");
 308     argList::validOptions.insert("scale", "scale factor");
 309
 310     argList args(argc, argv);
 311
 312     if (!args.check())
 313     {
 314         FatalError.exit();
 315     }
 316
 317     scalar scaleFactor = 1.0;
 318     args.optionReadIfPresent("scale", scaleFactor);
 319
 320 #   include "createTime.H"
 321
 322     fileName ansysFile(args.additionalArgs()[0]);
 323     ifstream ansysStream(ansysFile.c_str());
 324
 325     if (!ansysStream)
 326     {
 327         FatalErrorIn("ansysToFoam::main(int argc, char *argv[])")
 328             << args.executable()
 329             << ": file " << ansysFile << " not found"
 330             << exit(FatalError);
 331     }
 332
 333     yyFlexLexer lexer(&ansysStream);
 334     while (lexer.yylex() != 0)
 335     {}
 336
 337     Info<< "Creating points" << endl;
 338
 339     pointField points(slPoints.size());
 340
 341     label i = 0;
 342     forAllConstIter(SLList<point>, slPoints, pointIter)
 343     {
 344         // Scale points for the given scale factor
 345         points[i++] = scaleFactor * pointIter();
 346     }
 347
 348
 349     labelList pointMap(maxNodei+1);
 350
 351     i = 0;
 352     forAllConstIter(SLList<label>, slPointMap, pointMapIter)
 353     {
 354         pointMap[pointMapIter()] = i++;
 355     }
 356
 357     Info<< "Creating cells" << endl;
 358
 359     labelList cellMap(maxCelli+1);
 360
 361     i = 0;
 362     forAllConstIter(SLList<label>, slCellMap, cellMapIter)
 363     {
 364         cellMap[cellMapIter()] = i++;
 365     }
 366
 367
 368     const cellModel& hex = *(cellModeller::lookup("hex"));
 369     const cellModel& prism = *(cellModeller::lookup("prism"));
 370     const cellModel& pyr = *(cellModeller::lookup("pyr"));
 371     const cellModel& tet = *(cellModeller::lookup("tet"));
 372
 373     labelList labelsHex(8);
 374     labelList labelsPrism(6);
 375     labelList labelsPyramid(5);
 376     labelList labelsTet(4);
 377
 378     cellShapeList cellShapes(slCellLabels.size());
 379     label nCells = 0;
 380
 381     forAllConstIter(SLPtrList<labelList>, slCellLabels, cellIter)
 382     {
 383         if      // Tetrahedron
 384         (
 385             cellIter()[2] == cellIter()[3]
 386          && cellIter()[4] == cellIter()[5]
 387          && cellIter()[5] == cellIter()[6]
 388          && cellIter()[6] == cellIter()[7]
 389         )
 390         {
 391             labelsTet[0] = pointMap[cellIter()[0] ];
 392             labelsTet[1] = pointMap[cellIter()[1] ];
 393             labelsTet[2] = pointMap[cellIter()[2] ];
 394             labelsTet[3] = pointMap[cellIter()[4] ];
 395
 396             cellShapes[nCells++] = cellShape(tet, labelsTet);
 397         }
 398
 399         else if // Square-based pyramid
 400         (
 401             cellIter()[4] == cellIter()[5]
 402          && cellIter()[5] == cellIter()[6]
 403          && cellIter()[6] == cellIter()[7]
 404         )
 405         {
 406             labelsPyramid[0] = pointMap[cellIter()[0] ];
 407             labelsPyramid[1] = pointMap[cellIter()[1] ];
 408             labelsPyramid[2] = pointMap[cellIter()[2] ];
 409             labelsPyramid[3] = pointMap[cellIter()[3] ];
 410             labelsPyramid[4] = pointMap[cellIter()[4] ];
 411
 412             cellShapes[nCells++] = cellShape(pyr, labelsPyramid);
 413         }
 414
 415         else if // Triangular prism
 416         (
 417             cellIter()[2] == cellIter()[3]
 418          && cellIter()[6] == cellIter()[7]
 419         )
 420         {
 421             labelsPrism[0] = pointMap[cellIter()[0] ];
 422             labelsPrism[1] = pointMap[cellIter()[1] ];
 423             labelsPrism[2] = pointMap[cellIter()[2] ];
 424             labelsPrism[3] = pointMap[cellIter()[4] ];
 425             labelsPrism[4] = pointMap[cellIter()[5] ];
 426             labelsPrism[5] = pointMap[cellIter()[6] ];
 427
 428             cellShapes[nCells++] = cellShape(prism, labelsPrism);
 429         }
 430
 431         else // Hex
 432         {
 433             labelsHex[0] = pointMap[cellIter()[0] ];
 434             labelsHex[1] = pointMap[cellIter()[1] ];
 435             labelsHex[2] = pointMap[cellIter()[2] ];
 436             labelsHex[3] = pointMap[cellIter()[3] ];
 437             labelsHex[4] = pointMap[cellIter()[4] ];
 438             labelsHex[5] = pointMap[cellIter()[5] ];
 439             labelsHex[6] = pointMap[cellIter()[6] ];
 440             labelsHex[7] = pointMap[cellIter()[7] ];
 441
 442             cellShapes[nCells++] = cellShape(hex, labelsHex);
 443         }
 444     }
 445
 446
 447     const word defaultFacesName = "defaultFaces";
 448     word defaultFacesType = emptyPolyPatch::typeName;
 449
 450     // Create dummy mesh just to find out what are internal/external
 451     // faces
 452     autoPtr<polyMesh> dummyMesh
 453     (
 454         new polyMesh
 455         (
 456             IOobject
 457             (
 458                 "dummyMesh",
 459                 runTime.constant(),
 460                 runTime
 461             ),
 462             xferCopy(points),
 463             cellShapes,
 464             faceListList(0),
 465             wordList(0),
 466             wordList(0),
 467             defaultFacesName,
 468             defaultFacesType,
 469             wordList(0)
 470         )
 471     );
 472
 473
 474     // Warning: tet face order has changed between version 1.9.6 and 2.0
 475     //
 476     label faceIndex[7][6] =
 477     {
 478         {-1, -1, -1, -1, -1, -1}, // 0
 479         {-1, -1, -1, -1, -1, -1}, // 1
 480         {-1, -1, -1, -1, -1, -1}, // 2
 481         {-1, -1, -1, -1, -1, -1}, // 3
 482         {-1,  2,  0,  3,  1, -1}, // tet (version 2.0)
 483         { 0,  4,  3, -1,  2,  1}, // prism
 484         { 4,  2,  1,  3,  0,  5}, // hex
 485     };
 486
 487     Info<< "Creating boundary patches" << endl;
 488
 489     faceListList boundary(slPatchCells.size());
 490     wordList patchNames(slPatchCells.size());
 491
 492     forAll(slPatchCells, patchI)
 493     {
 494         SLList<face> patchFaces;
 495
 496         SLList<label>::iterator cellIter(slPatchCells[patchI].begin());
 497         SLList<label>::iterator faceIter(slPatchCellFaces[patchI].begin());
 498
 499         for
 500         (
 501             ;
 502             cellIter != slPatchCells[patchI].end()
 503          && faceIter != slPatchCellFaces[patchI].end();
 504             ++cellIter, ++faceIter
 505         )
 506         {
 507             const cellShape& shape = cellShapes[cellMap[cellIter()]];
 508
 509             patchFaces.append
 510             (
 511                 shape.faces()
 512                 [
 513                     faceIndex
 514                         [shape.nFaces()]
 515                         [faceIter()-1]
 516                 ]
 517             );
 518         }
 519
 520         boundary[patchI] = patchFaces;
 521         patchNames[patchI] = word("patch") + name(patchI + 1);
 522     }
 523
 524
 525     //
 526     // Lookup the face labels for all the boundary faces
 527     //
 528     labelListList boundaryFaceLabels(boundary.size());
 529     forAll(boundary, patchI)
 530     {
 531         const faceList& bFaces = boundary[patchI];
 532         labelList& bFaceLabels = boundaryFaceLabels[patchI];
 533         bFaceLabels.setSize(bFaces.size());
 534         forAll(bFaces, i)
 535         {
 536             bFaceLabels[i] = findFace(dummyMesh(), bFaces[i]);
 537         }
 538     }
 539
 540
 541     // Now split the boundary faces into external and internal faces. All
 542     // faces go into faceZones and external faces go into patches.
 543     List<faceList> patchFaces(slPatchCells.size());
 544     labelList patchNFaces(slPatchCells.size(), 0);
 545     forAll(boundary, patchI)
 546     {
 547         const faceList& bFaces = boundary[patchI];
 548         const labelList& bFaceLabels = boundaryFaceLabels[patchI];
 549
 550         patchFaces[patchI].setSize(bFaces.size());
 551
 552         forAll(bFaces, i)
 553         {
 554             if (!dummyMesh().isInternalFace(bFaceLabels[i]))
 555             {
 556                 patchFaces[patchI][patchNFaces[patchI]++] = bFaces[i];
 557             }
 558         }
 559         patchFaces[patchI].setSize(patchNFaces[patchI]);
 560
 561         Info<< "Patch " << patchI << " named " << patchNames[patchI]
 562             << ": " << boundary[patchI].size() << " faces" << endl;
 563     }
 564
 565     // We no longer need the dummyMesh
 566     dummyMesh.clear();
 567
 568
 569     Info<< "ansysToFoam: " << endl
 570         << "Ansys file format does not provide information about the type of "
 571         << "the patch (eg. wall, symmetry plane, cyclic etc)." << endl
 572         << "All the patches have been created "
 573         << "as type patch. Please reset after mesh conversion as necessary."
 574         << endl;
 575
 576     PtrList<dictionary> patchDicts;
 577
 578     preservePatchTypes
 579     (
 580         runTime,
 581         runTime.constant(),
 582         polyMesh::meshSubDir,
 583         patchNames,
 584         patchDicts,
 585         defaultFacesName,
 586         defaultFacesType
 587     );
 588
 589     // Add information to dictionary
 590     forAll(patchNames, patchI)
 591     {
 592         if (!patchDicts.set(patchI))
 593         {
 594             patchDicts.set(patchI, new dictionary());
 595         }
 596         // Add but not overwrite
 597         patchDicts[patchI].add("type", polyPatch::typeName, false);
 598     }
 599
 600
 601     polyMesh pShapeMesh
 602     (
 603         IOobject
 604         (
 605             polyMesh::defaultRegion,
 606             runTime.constant(),
 607             runTime
 608         ),
 609         xferMove(points),
 610         cellShapes,
 611         patchFaces,
 612         patchNames,
 613         patchDicts,
 614         defaultFacesName,
 615         defaultFacesType
 616     );
 617
 618
 619     if (cellTypes.size() > 0 || patchNames.size() > 0)
 620     {
 621         DynamicList<pointZone*> pz;
 622         DynamicList<faceZone*> fz;
 623         DynamicList<cellZone*> cz;
 624
 625         // FaceZones
 626         forAll(boundaryFaceLabels, patchI)
 627         {
 628             if (boundaryFaceLabels[patchI].size())
 629             {
 630                 // Re-do the boundaryFaceLabels since the boundary face
 631                 // labels will be different on the pShapeMesh.
 632                 const faceList& bFaces = boundary[patchI];
 633                 labelList& bFaceLabels = boundaryFaceLabels[patchI];
 634                 forAll(bFaceLabels, i)
 635                 {
 636                     bFaceLabels[i] = findFace(pShapeMesh, bFaces[i]);
 637                 }
 638
 639                 Info<< "Creating faceZone " <<  patchNames[patchI]
 640                     << " with " << bFaceLabels.size() << " faces" << endl;
 641
 642                 fz.append
 643                 (
 644                     new faceZone
 645                     (
 646                         patchNames[patchI],
 647                         bFaceLabels,
 648                         boolList(bFaceLabels.size(), false),
 649                         fz.size(),
 650                         pShapeMesh.faceZones()
 651                     )
 652                 );
 653             }
 654         }
 655
 656
 657         // CellZones
 658         labelList types = cellTypes.toc();
 659
 660         forAll(types, j)
 661         {
 662             label cellType = types[j];
 663
 664             // Pick up cells in zone
 665             DynamicList<label> addr;
 666
 667             SLList<label>::iterator cellMapIter = slCellMap.begin();
 668             SLList<label>::iterator typeIter = slCellType.begin();
 669
 670             for
 671             (
 672                 ;
 673                 typeIter != slCellType.end();
 674                 ++typeIter, ++cellMapIter
 675             )
 676             {
 677                 if (typeIter() == cellType)
 678                 {
 679                     addr.append(cellMap[cellMapIter()]);
 680                 }
 681             }
 682
 683             Info<< "Creating cellZone " << cellTypes[cellType]
 684                 << " with " << addr.size() << " cells" << endl;
 685
 686             cz.append
 687             (
 688                 new cellZone
 689                 (
 690                     cellTypes[cellType],
 691                     addr,
 692                     j,
 693                     pShapeMesh.cellZones()
 694                 )
 695             );
 696         }
 697
 698         pShapeMesh.addZones(pz, fz, cz);
 699     }
 700
 701
 702     // Set the precision of the points data to 10
 703     IOstream::defaultPrecision(max(10u, IOstream::defaultPrecision()));
 704
 705     Info<< "Writing polyMesh" << endl;
 706     pShapeMesh.write();
 707
 708     Info<< nl << "end" << endl;
 709     return 0;
 710 }
 711
 712
 713  /* ------------------------------------------------------------------------- *\
 714     ------ End of ansysToFoam.L
 715  \* ------------------------------------------------------------------------- */