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

   1 /*--------------------------------*- C++ -*----------------------------------*\
   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 Application
  25     ansysToFoam
  26
  27 Description
  28     Converts an ANSYS input mesh file, exported from I-DEAS,
  29     to OpenFOAM format.
  30
  31 \*---------------------------------------------------------------------------*/
  32
  33 %{
  34
  35 #undef yyFlexLexer
  36
  37  /* ------------------------------------------------------------------------- *\
  38    ------ local definitions
  39  \* ------------------------------------------------------------------------- */
  40
  41 #include <sstream>
  42 // For EOF only
  43 #include <cstdio>
  44
  45 #include "scalar.H"
  46 #include "IStringStream.H"
  47
  48 using namespace Foam;
  49
  50 #include "argList.H"
  51 #include "Time.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_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
 274 label findFace(const polyMesh& mesh, const face& f)
 275 {
 276     const labelList& pFaces = mesh.pointFaces()[f[0]];
 277
 278     forAll(pFaces, i)
 279     {
 280         label faceI = pFaces[i];
 281
 282         if (mesh.faces()[faceI] == f)
 283         {
 284             return faceI;
 285         }
 286     }
 287
 288     FatalErrorIn("findFace(const polyMesh&, const face&)")
 289         << "Cannot find a face matching " << f
 290         << exit(FatalError);
 291
 292     return -1;
 293 }
 294
 295
 296 int main(int argc, char *argv[])
 297 {
 298     argList::noParallel();
 299     argList::validArgs.append("ANSYS input file");
 300     argList::addOption
 301     (
 302         "scale",
 303         "factor",
 304         "geometry scaling factor - default is 1"
 305     );
 306
 307     argList args(argc, argv);
 308
 309     if (!args.check())
 310     {
 311         FatalError.exit();
 312     }
 313
 314     const scalar scaleFactor = args.optionLookupOrDefault("scale", 1.0);
 315
 316 #   include "createTime.H"
 317
 318     fileName ansysFile(args.additionalArgs()[0]);
 319     ifstream ansysStream(ansysFile.c_str());
 320
 321     if (!ansysStream)
 322     {
 323         FatalErrorIn("ansysToFoam::main(int argc, char *argv[])")
 324             << args.executable()
 325             << ": file " << ansysFile << " not found"
 326             << exit(FatalError);
 327     }
 328
 329     yyFlexLexer lexer(&ansysStream);
 330     while (lexer.yylex() != 0)
 331     {}
 332
 333     Info<< "Creating points" << endl;
 334
 335     pointField points(slPoints.size());
 336
 337     label i = 0;
 338     forAllConstIter(SLList<point>, slPoints, pointIter)
 339     {
 340         // Scale points for the given scale factor
 341         points[i++] = scaleFactor * pointIter();
 342     }
 343
 344
 345     labelList pointMap(maxNodei+1);
 346
 347     i = 0;
 348     forAllConstIter(SLList<label>, slPointMap, pointMapIter)
 349     {
 350         pointMap[pointMapIter()] = i++;
 351     }
 352
 353     Info<< "Creating cells" << endl;
 354
 355     labelList cellMap(maxCelli+1);
 356
 357     i = 0;
 358     forAllConstIter(SLList<label>, slCellMap, cellMapIter)
 359     {
 360         cellMap[cellMapIter()] = i++;
 361     }
 362
 363
 364     const cellModel& hex = *(cellModeller::lookup("hex"));
 365     const cellModel& prism = *(cellModeller::lookup("prism"));
 366     const cellModel& pyr = *(cellModeller::lookup("pyr"));
 367     const cellModel& tet = *(cellModeller::lookup("tet"));
 368
 369     labelList labelsHex(8);
 370     labelList labelsPrism(6);
 371     labelList labelsPyramid(5);
 372     labelList labelsTet(4);
 373
 374     cellShapeList cellShapes(slCellLabels.size());
 375     label nCells = 0;
 376
 377     forAllConstIter(SLPtrList<labelList>, slCellLabels, cellIter)
 378     {
 379         if      // Tetrahedron
 380         (
 381             cellIter()[2] == cellIter()[3]
 382          && cellIter()[4] == cellIter()[5]
 383          && cellIter()[5] == cellIter()[6]
 384          && cellIter()[6] == cellIter()[7]
 385         )
 386         {
 387             labelsTet[0] = pointMap[cellIter()[0] ];
 388             labelsTet[1] = pointMap[cellIter()[1] ];
 389             labelsTet[2] = pointMap[cellIter()[2] ];
 390             labelsTet[3] = pointMap[cellIter()[4] ];
 391
 392             cellShapes[nCells++] = cellShape(tet, labelsTet);
 393         }
 394
 395         else if // Square-based pyramid
 396         (
 397             cellIter()[4] == cellIter()[5]
 398          && cellIter()[5] == cellIter()[6]
 399          && cellIter()[6] == cellIter()[7]
 400         )
 401         {
 402             labelsPyramid[0] = pointMap[cellIter()[0] ];
 403             labelsPyramid[1] = pointMap[cellIter()[1] ];
 404             labelsPyramid[2] = pointMap[cellIter()[2] ];
 405             labelsPyramid[3] = pointMap[cellIter()[3] ];
 406             labelsPyramid[4] = pointMap[cellIter()[4] ];
 407
 408             cellShapes[nCells++] = cellShape(pyr, labelsPyramid);
 409         }
 410
 411         else if // Triangular prism
 412         (
 413             cellIter()[2] == cellIter()[3]
 414          && cellIter()[6] == cellIter()[7]
 415         )
 416         {
 417             labelsPrism[0] = pointMap[cellIter()[0] ];
 418             labelsPrism[1] = pointMap[cellIter()[1] ];
 419             labelsPrism[2] = pointMap[cellIter()[2] ];
 420             labelsPrism[3] = pointMap[cellIter()[4] ];
 421             labelsPrism[4] = pointMap[cellIter()[5] ];
 422             labelsPrism[5] = pointMap[cellIter()[6] ];
 423
 424             cellShapes[nCells++] = cellShape(prism, labelsPrism);
 425         }
 426
 427         else // Hex
 428         {
 429             labelsHex[0] = pointMap[cellIter()[0] ];
 430             labelsHex[1] = pointMap[cellIter()[1] ];
 431             labelsHex[2] = pointMap[cellIter()[2] ];
 432             labelsHex[3] = pointMap[cellIter()[3] ];
 433             labelsHex[4] = pointMap[cellIter()[4] ];
 434             labelsHex[5] = pointMap[cellIter()[5] ];
 435             labelsHex[6] = pointMap[cellIter()[6] ];
 436             labelsHex[7] = pointMap[cellIter()[7] ];
 437
 438             cellShapes[nCells++] = cellShape(hex, labelsHex);
 439         }
 440     }
 441
 442
 443     const word defaultFacesName = "defaultFaces";
 444     word defaultFacesType = emptyPolyPatch::typeName;
 445
 446     // Create dummy mesh just to find out what are internal/external
 447     // faces
 448     autoPtr<polyMesh> dummyMesh
 449     (
 450         new polyMesh
 451         (
 452             IOobject
 453             (
 454                 "dummyMesh",
 455                 runTime.constant(),
 456                 runTime
 457             ),
 458             xferCopy(points),
 459             cellShapes,
 460             faceListList(0),
 461             wordList(0),
 462             wordList(0),
 463             defaultFacesName,
 464             defaultFacesType,
 465             wordList(0)
 466         )
 467     );
 468
 469
 470     // Warning: tet face order has changed between version 1.9.6 and 2.0
 471     //
 472     label faceIndex[7][6] =
 473     {
 474         {-1, -1, -1, -1, -1, -1}, // 0
 475         {-1, -1, -1, -1, -1, -1}, // 1
 476         {-1, -1, -1, -1, -1, -1}, // 2
 477         {-1, -1, -1, -1, -1, -1}, // 3
 478         { 3,  2,  0, -1,  1, -1}, // tet (version 2.0)
 479         { 0,  4,  3, -1,  2,  1}, // prism
 480         { 4,  2,  1,  3,  0,  5}, // hex
 481     };
 482
 483     Info<< "Creating boundary patches" << endl;
 484
 485     faceListList boundary(slPatchCells.size());
 486     wordList patchNames(slPatchCells.size());
 487
 488     forAll(slPatchCells, patchI)
 489     {
 490         SLList<face> patchFaces;
 491
 492         SLList<label>::iterator cellIter(slPatchCells[patchI].begin());
 493         SLList<label>::iterator faceIter(slPatchCellFaces[patchI].begin());
 494
 495         for
 496         (
 497             ;
 498             cellIter != slPatchCells[patchI].end()
 499          && faceIter != slPatchCellFaces[patchI].end();
 500             ++cellIter, ++faceIter
 501         )
 502         {
 503             const cellShape& shape = cellShapes[cellMap[cellIter()]];
 504
 505             patchFaces.append
 506             (
 507                 shape.faces()
 508                 [
 509                     faceIndex
 510                         [shape.nFaces()]
 511                         [faceIter()-1]
 512                 ]
 513             );
 514         }
 515
 516         boundary[patchI] = patchFaces;
 517         patchNames[patchI] = word("patch") + name(patchI + 1);
 518     }
 519
 520
 521     //
 522     // Lookup the face labels for all the boundary faces
 523     //
 524     labelListList boundaryFaceLabels(boundary.size());
 525     forAll(boundary, patchI)
 526     {
 527         const faceList& bFaces = boundary[patchI];
 528         labelList& bFaceLabels = boundaryFaceLabels[patchI];
 529         bFaceLabels.setSize(bFaces.size());
 530         forAll(bFaces, i)
 531         {
 532             bFaceLabels[i] = findFace(dummyMesh(), bFaces[i]);
 533         }
 534     }
 535
 536
 537     // Now split the boundary faces into external and internal faces. All
 538     // faces go into faceZones and external faces go into patches.
 539     List<faceList> patchFaces(slPatchCells.size());
 540     labelList patchNFaces(slPatchCells.size(), 0);
 541     forAll(boundary, patchI)
 542     {
 543         const faceList& bFaces = boundary[patchI];
 544         const labelList& bFaceLabels = boundaryFaceLabels[patchI];
 545
 546         patchFaces[patchI].setSize(bFaces.size());
 547
 548         forAll(bFaces, i)
 549         {
 550             if (!dummyMesh().isInternalFace(bFaceLabels[i]))
 551             {
 552                 patchFaces[patchI][patchNFaces[patchI]++] = bFaces[i];
 553             }
 554         }
 555         patchFaces[patchI].setSize(patchNFaces[patchI]);
 556
 557         Info<< "Patch " << patchI << " named " << patchNames[patchI]
 558             << ": " << boundary[patchI].size() << " faces" << endl;
 559     }
 560
 561     // We no longer need the dummyMesh
 562     dummyMesh.clear();
 563
 564
 565     Info<< "ansysToFoam: " << endl
 566         << "Ansys file format does not provide information about the type of "
 567         << "the patch (eg. wall, symmetry plane, cyclic etc)." << endl
 568         << "All the patches have been created "
 569         << "as type patch. Please reset after mesh conversion as necessary."
 570         << endl;
 571
 572     PtrList<dictionary> patchDicts;
 573
 574     preservePatchTypes
 575     (
 576         runTime,
 577         runTime.constant(),
 578         polyMesh::meshSubDir,
 579         patchNames,
 580         patchDicts,
 581         defaultFacesName,
 582         defaultFacesType
 583     );
 584
 585     // Add information to dictionary
 586     forAll(patchNames, patchI)
 587     {
 588         if (!patchDicts.set(patchI))
 589         {
 590             patchDicts.set(patchI, new dictionary());
 591         }
 592         // Add but not overwrite
 593         patchDicts[patchI].add("type", polyPatch::typeName, false);
 594     }
 595
 596
 597     polyMesh pShapeMesh
 598     (
 599         IOobject
 600         (
 601             polyMesh::defaultRegion,
 602             runTime.constant(),
 603             runTime
 604         ),
 605         xferMove(points),
 606         cellShapes,
 607         patchFaces,
 608         patchNames,
 609         patchDicts,
 610         defaultFacesName,
 611         defaultFacesType
 612     );
 613
 614
 615     if (cellTypes.size() > 0 || patchNames.size() > 0)
 616     {
 617         DynamicList<pointZone*> pz;
 618         DynamicList<faceZone*> fz;
 619         DynamicList<cellZone*> cz;
 620
 621         // FaceZones
 622         forAll(boundaryFaceLabels, patchI)
 623         {
 624             if (boundaryFaceLabels[patchI].size())
 625             {
 626                 // Re-do the boundaryFaceLabels since the boundary face
 627                 // labels will be different on the pShapeMesh.
 628                 const faceList& bFaces = boundary[patchI];
 629                 labelList& bFaceLabels = boundaryFaceLabels[patchI];
 630                 forAll(bFaceLabels, i)
 631                 {
 632                     bFaceLabels[i] = findFace(pShapeMesh, bFaces[i]);
 633                 }
 634
 635                 Info<< "Creating faceZone " <<  patchNames[patchI]
 636                     << " with " << bFaceLabels.size() << " faces" << endl;
 637
 638                 fz.append
 639                 (
 640                     new faceZone
 641                     (
 642                         patchNames[patchI],
 643                         bFaceLabels,
 644                         boolList(bFaceLabels.size(), false),
 645                         fz.size(),
 646                         pShapeMesh.faceZones()
 647                     )
 648                 );
 649             }
 650         }
 651
 652
 653         // CellZones
 654         labelList types = cellTypes.sortedToc();
 655
 656         forAll(types, j)
 657         {
 658             label cellType = types[j];
 659
 660             // Pick up cells in zone
 661             DynamicList<label> addr;
 662
 663             SLList<label>::iterator cellMapIter = slCellMap.begin();
 664             SLList<label>::iterator typeIter = slCellType.begin();
 665
 666             for
 667             (
 668                 ;
 669                 typeIter != slCellType.end();
 670                 ++typeIter, ++cellMapIter
 671             )
 672             {
 673                 if (typeIter() == cellType)
 674                 {
 675                     addr.append(cellMap[cellMapIter()]);
 676                 }
 677             }
 678
 679             Info<< "Creating cellZone " << cellTypes[cellType]
 680                 << " with " << addr.size() << " cells" << endl;
 681
 682             cz.append
 683             (
 684                 new cellZone
 685                 (
 686                     cellTypes[cellType],
 687                     addr,
 688                     j,
 689                     pShapeMesh.cellZones()
 690                 )
 691             );
 692         }
 693
 694         pShapeMesh.addZones(pz, fz, cz);
 695     }
 696
 697
 698     // Set the precision of the points data to 10
 699     IOstream::defaultPrecision(10);
 700
 701     Info<< "Writing polyMesh" << endl;
 702     pShapeMesh.write();
 703
 704     Info<< nl << "end" << endl;
 705     return 0;
 706 }
 707
 708
 709  /* ------------------------------------------------------------------------- *\
 710     ------ End of ansysToFoam.L
 711  \* ------------------------------------------------------------------------- */