applications/utilities/mesh/conversion/plot3dToFoam/plot3dToFoam.C

   1 /*---------------------------------------------------------------------------*\
   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 Description
  25     Plot3d mesh (ascii/formatted format) converter.
  26
  27     Work in progress! Handles ascii multiblock (and optionally singleBlock)
  28     format.
  29     By default expects blanking. Use -noBlank if none.
  30     Use -2D @a thickness if 2D.
  31
  32     Niklas Nordin has experienced a problem with lefthandedness of the blocks.
  33     The code should detect this automatically - see hexBlock::readPoints but
  34     if this goes wrong just set the blockHandedness_ variable to 'right'
  35     always.
  36
  37 \*---------------------------------------------------------------------------*/
  38
  39 #include "argList.H"
  40 #include "objectRegistry.H"
  41 #include "foamTime.H"
  42 #include "IFstream.H"
  43 #include "hexBlock.H"
  44 #include "polyMesh.H"
  45 #include "wallPolyPatch.H"
  46 #include "symmetryPolyPatch.H"
  47 #include "preservePatchTypes.H"
  48 #include "cellShape.H"
  49 #include "cellModeller.H"
  50 #include "mergePoints.H"
  51
  52 using namespace Foam;
  53
  54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  55 // Main program:
  56
  57 int main(int argc, char *argv[])
  58 {
  59     argList::noParallel();
  60     argList::validArgs.append("PLOT3D geom file");
  61     argList::validOptions.insert("scale", "scale factor");
  62     argList::validOptions.insert("noBlank", "");
  63     argList::validOptions.insert("singleBlock", "");
  64     argList::validOptions.insert("2D", "thickness");
  65
  66     argList args(argc, argv);
  67
  68     if (!args.check())
  69     {
  70          FatalError.exit();
  71     }
  72
  73     scalar scaleFactor = 1.0;
  74     args.optionReadIfPresent("scale", scaleFactor);
  75
  76     bool readBlank = !args.optionFound("noBlank");
  77     bool singleBlock = args.optionFound("singleBlock");
  78     scalar twoDThickness = -1;
  79     if (args.optionReadIfPresent("2D", twoDThickness))
  80     {
  81         Info<< "Reading 2D case by extruding points by " << twoDThickness
  82             << " in z direction." << nl << endl;
  83     }
  84
  85
  86 #   include "createTime.H"
  87
  88     IFstream plot3dFile(args.additionalArgs()[0]);
  89
  90     // Read the plot3d information using a fixed format reader.
  91     // Comments in the file are in C++ style, so the stream parser will remove
  92     // them with no intervention
  93     label nblock;
  94
  95     if (singleBlock)
  96     {
  97         nblock = 1;
  98     }
  99     else
 100     {
 101         plot3dFile >> nblock;
 102     }
 103
 104     Info<< "Reading " << nblock << " blocks" << endl;
 105
 106     PtrList<hexBlock> blocks(nblock);
 107
 108     {
 109         label nx, ny, nz;
 110
 111         forAll (blocks, blockI)
 112         {
 113             if (twoDThickness > 0)
 114             {
 115                 // Fake second set of points (done in readPoints below)
 116                 plot3dFile >> nx >> ny;
 117                 nz = 2;
 118             }
 119             else
 120             {
 121                 plot3dFile >> nx >> ny >> nz;
 122             }
 123
 124             Info<< "block " << blockI << " nx:" << nx
 125                 << " ny:" << ny << " nz:" << nz << endl;
 126
 127             blocks.set(blockI, new hexBlock(nx, ny, nz));
 128         }
 129     }
 130
 131     Info<< "Reading block points" << endl;
 132     label sumPoints(0);
 133     label nMeshCells(0);
 134
 135     forAll (blocks, blockI)
 136     {
 137         Info<< "block " << blockI << ":" << nl;
 138         blocks[blockI].readPoints(readBlank, twoDThickness, plot3dFile);
 139         sumPoints += blocks[blockI].nBlockPoints();
 140         nMeshCells += blocks[blockI].nBlockCells();
 141         Info<< nl;
 142     }
 143
 144     pointField points(sumPoints);
 145     labelList blockOffsets(blocks.size());
 146     sumPoints = 0;
 147     forAll (blocks, blockI)
 148     {
 149         const pointField& blockPoints = blocks[blockI].points();
 150         blockOffsets[blockI] = sumPoints;
 151         forAll (blockPoints, i)
 152         {
 153             points[sumPoints++] = blockPoints[i];
 154         }
 155     }
 156
 157     // From old to new master point
 158     labelList oldToNew;
 159     pointField newPoints;
 160
 161     // Merge points
 162     mergePoints
 163     (
 164         points,
 165         SMALL,
 166         false,
 167         oldToNew,
 168         newPoints
 169     );
 170
 171     Info<< "Merged points within " << SMALL << " distance. Merged from "
 172         << oldToNew.size() << " down to " << newPoints.size()
 173         << " points." << endl;
 174
 175     // Scale the points
 176     if (scaleFactor > 1.0 + SMALL || scaleFactor < 1.0 - SMALL)
 177     {
 178         newPoints *= scaleFactor;
 179     }
 180
 181     Info<< "Creating cells" << endl;
 182
 183     cellShapeList cellShapes(nMeshCells);
 184
 185     const cellModel& hex = *(cellModeller::lookup("hex"));
 186
 187     label nCreatedCells = 0;
 188
 189     forAll (blocks, blockI)
 190     {
 191         labelListList curBlockCells = blocks[blockI].blockCells();
 192
 193         forAll (curBlockCells, blockCellI)
 194         {
 195             labelList cellPoints(curBlockCells[blockCellI].size());
 196
 197             forAll (cellPoints, pointI)
 198             {
 199                 cellPoints[pointI] =
 200                     oldToNew
 201                     [
 202                         curBlockCells[blockCellI][pointI]
 203                       + blockOffsets[blockI]
 204                     ];
 205             }
 206
 207             // Do automatic collapse from hex.
 208             cellShapes[nCreatedCells] = cellShape(hex, cellPoints, true);
 209
 210             nCreatedCells++;
 211         }
 212     }
 213
 214     Info<< "Creating boundary patches" << endl;
 215
 216     faceListList boundary(0);
 217     wordList patchNames(0);
 218     wordList patchTypes(0);
 219     word defaultFacesName = "defaultFaces";
 220     word defaultFacesType = wallPolyPatch::typeName;
 221     wordList patchPhysicalTypes(0);
 222
 223     polyMesh pShapeMesh
 224     (
 225         IOobject
 226         (
 227             polyMesh::defaultRegion,
 228             runTime.constant(),
 229             runTime
 230         ),
 231         xferMove(newPoints),
 232         cellShapes,
 233         boundary,
 234         patchNames,
 235         patchTypes,
 236         defaultFacesName,
 237         defaultFacesType,
 238         patchPhysicalTypes
 239     );
 240
 241     // Set the precision of the points data to 10
 242     IOstream::defaultPrecision(10);
 243
 244     Info<< "Writing polyMesh" << endl;
 245     pShapeMesh.write();
 246
 247     Info<< "End\n" << endl;
 248
 249     return 0;
 250 }
 251
 252
 253 // ************************************************************************* //