FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

STARCDMeshReader.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "STARCDMeshReader.H"
00027 #include <OpenFOAM/cyclicPolyPatch.H>
00028 #include <OpenFOAM/emptyPolyPatch.H>
00029 #include <OpenFOAM/wallPolyPatch.H>
00030 #include <OpenFOAM/symmetryPolyPatch.H>
00031 #include <OpenFOAM/cellModeller.H>
00032 #include <OpenFOAM/ListOps.H>
00033 #include <OpenFOAM/IFstream.H>
00034 #include <OpenFOAM/IOMap.H>
00035 
00036 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00037 
00038 const char* Foam::meshReaders::STARCD::defaultBoundaryName =
00039     "Default_Boundary_Region";
00040 
00041 const char* Foam::meshReaders::STARCD::defaultSolidBoundaryName =
00042     "Default_Boundary_Solid";
00043 
00044 bool Foam::meshReaders::STARCD::keepSolids = false;
00045 
00046 const int Foam::meshReaders::STARCD::starToFoamFaceAddr[4][6] =
00047 {
00048     { 4, 5, 2, 3, 0, 1 },     // 11 = pro-STAR hex
00049     { 0, 1, 4, -1, 2, 3 },    // 12 = pro-STAR prism
00050     { 3, -1, 2, -1, 1, 0 },   // 13 = pro-STAR tetra
00051     { 0, -1, 4, 2, 1, 3 }     // 14 = pro-STAR pyramid
00052 };
00053 
00054 
00055 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00056 
00057 void Foam::meshReaders::STARCD::readToNewline(IFstream& is)
00058 {
00059     char ch = '\n';
00060     do
00061     {
00062         (is).get(ch);
00063     }
00064     while ((is) && ch != '\n');
00065 }
00066 
00067 
00068 bool Foam::meshReaders::STARCD::readHeader(IFstream& is, word fileSignature)
00069 {
00070     if (!is.good())
00071     {
00072         FatalErrorIn("meshReaders::STARCD::readHeader()")
00073             << "cannot read " << fileSignature  << "  " << is.name()
00074             << abort(FatalError);
00075     }
00076 
00077     word header;
00078     label majorVersion;
00079 
00080     is >> header;
00081     is >> majorVersion;
00082 
00083     // skip the rest of the line
00084     readToNewline(is);
00085 
00086     // add other checks ...
00087     if (header != fileSignature)
00088     {
00089         Info<< "header mismatch " << fileSignature << "  " << is.name()
00090             << endl;
00091     }
00092 
00093     return true;
00094 }
00095 
00096 
00097 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00098 
00099 void Foam::meshReaders::STARCD::readAux(const objectRegistry& registry)
00100 {
00101     boundaryRegion_.readDict(registry);
00102     cellTable_.readDict(registry);
00103 }
00104 
00105 
00106 // read in the points from the .vrt file
00107 //
00108 /*---------------------------------------------------------------------------*\
00109 Line 1:
00110   PROSTAR_VERTEX [newline]
00111 
00112 Line 2:
00113   <version> 0 0 0 0 0 0 0 [newline]
00114 
00115 Body:
00116   <vertexId>  <x>  <y>  <z> [newline]
00117 
00118 \*---------------------------------------------------------------------------*/
00119 void Foam::meshReaders::STARCD::readPoints
00120 (
00121     const fileName& inputName,
00122     const scalar scaleFactor
00123 )
00124 {
00125     const word fileSignature = "PROSTAR_VERTEX";
00126     label nPoints = 0, maxId = 0;
00127 
00128     // Pass 1:
00129     // get # points and maximum vertex label
00130     {
00131         IFstream is(inputName);
00132         readHeader(is, fileSignature);
00133 
00134         label lineLabel;
00135         scalar x, y, z;
00136 
00137         while ((is >> lineLabel).good())
00138         {
00139             nPoints++;
00140             maxId = max(maxId, lineLabel);
00141             is >> x >> y >> z;
00142         }
00143     }
00144 
00145     Info<< "Number of points  = " << nPoints << endl;
00146 
00147     // set sizes and reset to invalid values
00148 
00149     points_.setSize(nPoints);
00150     mapToFoamPointId_.setSize(maxId+1);
00151 
00152     //- original Point number for a given vertex
00153     // might need again in the future
00156 
00157     mapToFoamPointId_ = -1;
00158 
00159     // Pass 2:
00160     // construct pointList and conversion table
00161     // from Star vertex numbers to Foam point labels
00162     if (nPoints > 0)
00163     {
00164         IFstream is(inputName);
00165         readHeader(is, fileSignature);
00166 
00167         label lineLabel;
00168 
00169         label pointI = 0;
00170         while ((is >> lineLabel).good())
00171         {
00172             is  >> points_[pointI].x()
00173                 >> points_[pointI].y()
00174                 >> points_[pointI].z();
00175 
00176             // might need again in the future
00178             mapToFoamPointId_[lineLabel] = pointI;
00179             pointI++;
00180         }
00181 
00182         if (nPoints > pointI)
00183         {
00184             nPoints = pointI;
00185             points_.setSize(nPoints);
00186             // might need again in the future
00188         }
00189 
00190         if (scaleFactor > 1.0 + SMALL || scaleFactor < 1.0 - SMALL)
00191         {
00192             points_ *= scaleFactor;
00193         }
00194     }
00195     else
00196     {
00197         FatalErrorIn("meshReaders::STARCD::readPoints()")
00198             << "no points in file " << inputName
00199             << abort(FatalError);
00200     }
00201 
00202 }
00203 
00204 
00205 // read in the cells from the .cel file
00206 //
00207 /*---------------------------------------------------------------------------*\
00208 Line 1:
00209   PROSTAR_CELL [newline]
00210 
00211 Line 2:
00212   <version> 0 0 0 0 0 0 0 [newline]
00213 
00214 Body:
00215   <cellId>  <shapeId>  <nLabels>  <cellTableId>  <typeId> [newline]
00216   <cellId>  <int1> .. <int8>
00217   <cellId>  <int9> .. <int16>
00218 
00219  with shapeId:
00220  *   1 = point
00221  *   2 = line
00222  *   3 = shell
00223  *  11 = hexa
00224  *  12 = prism
00225  *  13 = tetra
00226  *  14 = pyramid
00227  * 255 = polyhedron
00228 
00229  with typeId
00230  *   1 = fluid
00231  *   2 = solid
00232  *   3 = baffle
00233  *   4 = shell
00234  *   5 = line
00235  *   6 = point
00236 
00237 For primitive cell shapes, the number of vertices will never exceed 8 (hexa)
00238 and corresponds to <nLabels>.
00239 For polyhedral, <nLabels> includess an index table comprising beg/end pairs
00240 for each cell face.
00241 
00242 Strictly speaking, we only need the cellModeller for adding boundaries.
00243 \*---------------------------------------------------------------------------*/
00244 
00245 void Foam::meshReaders::STARCD::readCells(const fileName& inputName)
00246 {
00247     const word fileSignature = "PROSTAR_CELL";
00248     label nFluids = 0, nSolids = 0, nBaffles = 0, nShells = 0;
00249     label maxId = 0;
00250 
00251     bool unknownVertices = false;
00252 
00253 
00254     // Pass 1:
00255     // count nFluids, nSolids, nBaffle, nShell and maxId
00256     // also see if polyhedral cells were used
00257     {
00258         IFstream is(inputName);
00259         readHeader(is, fileSignature);
00260 
00261         label lineLabel, shapeId, nLabels, cellTableId, typeId;
00262 
00263         while ((is >> lineLabel).good())
00264         {
00265             label starCellId = lineLabel;
00266             is  >> shapeId
00267                 >> nLabels
00268                 >> cellTableId
00269                 >> typeId;
00270 
00271             // skip the rest of the line
00272             readToNewline(is);
00273 
00274             // max 8 indices per line
00275             while (nLabels > 0)
00276             {
00277                 readToNewline(is);
00278                 nLabels -= 8;
00279             }
00280 
00281             if (typeId == starcdFluidType)
00282             {
00283                 nFluids++;
00284                 maxId = max(maxId, starCellId);
00285 
00286                 if (!cellTable_.found(cellTableId))
00287                 {
00288                     cellTable_.setName(cellTableId);
00289                     cellTable_.setMaterial(cellTableId, "fluid");
00290                 }
00291             }
00292             else if (typeId == starcdSolidType)
00293             {
00294                 nSolids++;
00295                 if (keepSolids)
00296                 {
00297                     maxId = max(maxId, starCellId);
00298                 }
00299 
00300                 if (!cellTable_.found(cellTableId))
00301                 {
00302                     cellTable_.setName(cellTableId);
00303                     cellTable_.setMaterial(cellTableId, "solid");
00304                 }
00305 
00306             }
00307             else if (typeId == starcdBaffleType)
00308             {
00309                 // baffles have no cellTable entry
00310                 nBaffles++;
00311                 maxId = max(maxId, starCellId);
00312             }
00313             else if (typeId == starcdShellType)
00314             {
00315                 nShells++;
00316                 if (!cellTable_.found(cellTableId))
00317                 {
00318                     cellTable_.setName(cellTableId);
00319                     cellTable_.setMaterial(cellTableId, "shell");
00320                 }
00321             }
00322 
00323         }
00324     }
00325 
00326     Info<< "Number of fluids  = " << nFluids << nl
00327         << "Number of baffles = " << nBaffles << nl;
00328     if (keepSolids)
00329     {
00330         Info<< "Number of solids  = " << nSolids << nl;
00331     }
00332     else
00333     {
00334         Info<< "Ignored   solids  = " << nSolids << nl;
00335     }
00336     Info<< "Ignored   shells  = " << nShells << endl;
00337 
00338 
00339     label nCells;
00340     if (keepSolids)
00341     {
00342         nCells = nFluids + nSolids;
00343     }
00344     else
00345     {
00346         nCells = nFluids;
00347     }
00348 
00349     cellFaces_.setSize(nCells);
00350     cellShapes_.setSize(nCells);
00351     cellTableId_.setSize(nCells);
00352 
00353     // information for the interfaces
00354     baffleFaces_.setSize(nBaffles);
00355 
00356     // extra space for baffles
00357     origCellId_.setSize(nCells + nBaffles);
00358     mapToFoamCellId_.setSize(maxId+1);
00359     mapToFoamCellId_ = -1;
00360 
00361 
00362     // avoid undefined shapes for polyhedra
00363     cellShape genericShape(*unknownModel, labelList(0));
00364 
00365     // Pass 2:
00366     // construct cellFaces_ and possibly cellShapes_
00367     if (nCells <= 0)
00368     {
00369         FatalErrorIn("meshReaders::STARCD::readCells()")
00370             << "no cells in file " << inputName
00371             << abort(FatalError);
00372     }
00373     else
00374     {
00375         IFstream is(inputName);
00376         readHeader(is, fileSignature);
00377 
00378         labelList starLabels(64);
00379         label lineLabel, shapeId, nLabels, cellTableId, typeId;
00380 
00381         label cellI = 0;
00382         label baffleI = 0;
00383 
00384         while ((is >> lineLabel).good())
00385         {
00386             label starCellId = lineLabel;
00387             is  >> shapeId
00388                 >> nLabels
00389                 >> cellTableId
00390                 >> typeId;
00391 
00392             if (nLabels > starLabels.size())
00393             {
00394                 starLabels.setSize(nLabels);
00395             }
00396             starLabels = -1;
00397 
00398             // read indices - max 8 per line
00399             for (label i = 0; i < nLabels; ++i)
00400             {
00401                 if ((i % 8) == 0)
00402                 {
00403                     is >> lineLabel;
00404                 }
00405                 is >> starLabels[i];
00406             }
00407 
00408             // skip solid cells
00409             if (typeId == starcdSolidType && !keepSolids)
00410             {
00411                 continue;
00412             }
00413 
00414             // determine the foam cell shape
00415             const cellModel* curModelPtr = NULL;
00416 
00417             // fluid/solid cells
00418             switch (shapeId)
00419             {
00420                 case starcdHex:
00421                     curModelPtr = hexModel;
00422                     break;
00423                 case starcdPrism:
00424                     curModelPtr = prismModel;
00425                     break;
00426                 case starcdTet:
00427                     curModelPtr = tetModel;
00428                     break;
00429                 case starcdPyr:
00430                     curModelPtr = pyrModel;
00431                     break;
00432             }
00433 
00434             if (curModelPtr)
00435             {
00436                 // primitive cell - use shapes
00437 
00438                 // convert orig vertex Id to point label
00439                 bool isBad = false;
00440                 for (label i=0; i < nLabels; ++i)
00441                 {
00442                     label pointId = mapToFoamPointId_[starLabels[i]];
00443                     if (pointId < 0)
00444                     {
00445                         Info<< "Cells inconsistent with vertex file. "
00446                             << "Star vertex " << starLabels[i]
00447                             << " does not exist" << endl;
00448                         isBad = true;
00449                         unknownVertices = true;
00450                     }
00451                     starLabels[i] = pointId;
00452                 }
00453 
00454                 if (isBad)
00455                 {
00456                     continue;
00457                 }
00458 
00459                 // record original cell number and lookup
00460                 origCellId_[cellI] = starCellId;
00461                 mapToFoamCellId_[starCellId] = cellI;
00462 
00463                 cellTableId_[cellI] = cellTableId;
00464                 cellShapes_[cellI] = cellShape
00465                 (
00466                     *curModelPtr,
00467                     SubList<label>(starLabels, nLabels)
00468                 );
00469 
00470                 cellFaces_[cellI] = cellShapes_[cellI].faces();
00471                 cellI++;
00472             }
00473             else if (shapeId == starcdPoly)
00474             {
00475                 // polyhedral cell
00476                 label nFaces = starLabels[0] - 1;
00477 
00478                 // convert orig vertex id to point label
00479                 // start with offset (skip the index table)
00480                 bool isBad = false;
00481                 for (label i=starLabels[0]; i < nLabels; ++i)
00482                 {
00483                     label pointId = mapToFoamPointId_[starLabels[i]];
00484                     if (pointId < 0)
00485                     {
00486                         Info<< "Cells inconsistent with vertex file. "
00487                             << "Star vertex " << starLabels[i]
00488                             << " does not exist" << endl;
00489                         isBad = true;
00490                         unknownVertices = true;
00491                     }
00492                     starLabels[i] = pointId;
00493                 }
00494 
00495                 if (isBad)
00496                 {
00497                     continue;
00498                 }
00499 
00500                 // traverse beg/end indices
00501                 faceList faces(nFaces);
00502                 label faceI = 0;
00503                 for (label i=0; i < nFaces; ++i)
00504                 {
00505                     label beg = starLabels[i];
00506                     label n   = starLabels[i+1] - beg;
00507 
00508                     face f
00509                     (
00510                         SubList<label>(starLabels, n, beg)
00511                     );
00512 
00513                     f.collapse();
00514 
00515                     // valid faces only
00516                     if (f.size() >= 3)
00517                     {
00518                         faces[faceI++] = f;
00519                     }
00520                 }
00521 
00522                 if (nFaces > faceI)
00523                 {
00524                     Info<< "star cell " << starCellId << " has "
00525                         << (nFaces - faceI)
00526                         << " empty faces - could cause boundary "
00527                         << "addressing problems"
00528                         << endl;
00529 
00530                     nFaces = faceI;
00531                     faces.setSize(nFaces);
00532                 }
00533 
00534                 if (nFaces < 4)
00535                 {
00536                     FatalErrorIn("meshReaders::STARCD::readCells()")
00537                         << "star cell " << starCellId << " has " << nFaces
00538                         << abort(FatalError);
00539                 }
00540 
00541                 // record original cell number and lookup
00542                 origCellId_[cellI] = starCellId;
00543                 mapToFoamCellId_[starCellId] = cellI;
00544 
00545                 cellTableId_[cellI] = cellTableId;
00546                 cellShapes_[cellI]  = genericShape;
00547                 cellFaces_[cellI]   = faces;
00548                 cellI++;
00549             }
00550             else if (typeId == starcdBaffleType)
00551             {
00552                 // baffles
00553 
00554                 // convert orig vertex id to point label
00555                 bool isBad = false;
00556                 for (label i=0; i < nLabels; ++i)
00557                 {
00558                     label pointId = mapToFoamPointId_[starLabels[i]];
00559                     if (pointId < 0)
00560                     {
00561                         Info<< "Baffles inconsistent with vertex file. "
00562                             << "Star vertex " << starLabels[i]
00563                             << " does not exist" << endl;
00564                         isBad = true;
00565                         unknownVertices = true;
00566                     }
00567                     starLabels[i] = pointId;
00568                 }
00569 
00570                 if (isBad)
00571                 {
00572                     continue;
00573                 }
00574 
00575 
00576                 face f
00577                 (
00578                     SubList<label>(starLabels, nLabels)
00579                 );
00580 
00581                 f.collapse();
00582 
00583                 // valid faces only
00584                 if (f.size() >= 3)
00585                 {
00586                     baffleFaces_[baffleI] = f;
00587                     // insert lookup addressing in normal list
00588                     mapToFoamCellId_[starCellId]  = nCells + baffleI;
00589                     origCellId_[nCells + baffleI] = starCellId;
00590                     baffleI++;
00591                 }
00592             }
00593         }
00594 
00595         baffleFaces_.setSize(baffleI);
00596     }
00597 
00598     if (unknownVertices)
00599     {
00600         FatalErrorIn("meshReaders::STARCD::readCells()")
00601             << "cells with unknown vertices"
00602             << abort(FatalError);
00603     }
00604 
00605     // truncate lists
00606 
00607 #ifdef DEBUG_READING
00608     Info<< "CELLS READ" << endl;
00609 #endif
00610 
00611     // cleanup
00612     mapToFoamPointId_.clear();
00613 }
00614 
00615 
00616 // read in the boundaries from the .bnd file
00617 //
00618 /*---------------------------------------------------------------------------*\
00619 Line 1:
00620   PROSTAR_BOUNDARY [newline]
00621 
00622 Line 2:
00623   <version> 0 0 0 0 0 0 0 [newline]
00624 
00625 Body:
00626   <boundId>  <cellId>  <cellFace>  <regionId>  0  <boundaryType> [newline]
00627 
00628 where boundaryType is truncated to 4 characters from one of the following:
00629 INLET
00630 PRESSSURE
00631 OUTLET
00632 BAFFLE
00633 etc,
00634 \*---------------------------------------------------------------------------*/
00635 
00636 void Foam::meshReaders::STARCD::readBoundary(const fileName& inputName)
00637 {
00638     const word fileSignature = "PROSTAR_BOUNDARY";
00639     label nPatches = 0, nFaces = 0, nBafflePatches = 0, maxId = 0;
00640     label lineLabel, starCellId, cellFaceId, starRegion, configNumber;
00641     word patchType;
00642 
00643     labelList mapToFoamPatchId(1000, -1);
00644     labelList nPatchFaces(1000, 0);
00645     labelList origRegion(1000, 0);
00646     patchTypes_.setSize(1000);
00647 
00648     // this is what we seem to need
00649     // these MUST correspond to starToFoamFaceAddr
00650     //
00651     Map<label> faceLookupIndex;
00652 
00653     faceLookupIndex.insert(hexModel->index(), 0);
00654     faceLookupIndex.insert(prismModel->index(), 1);
00655     faceLookupIndex.insert(tetModel->index(), 2);
00656     faceLookupIndex.insert(pyrModel->index(), 3);
00657 
00658     // Pass 1:
00659     // collect
00660     // no. of faces (nFaces), no. of patches (nPatches)
00661     // and for each of these patches the number of faces
00662     // (nPatchFaces[patchLabel])
00663     //
00664     // and a conversion table from Star regions to (Foam) patchLabels
00665     //
00666     // additionally note the no. of baffle patches (nBafflePatches)
00667     // so that we sort these to the end of the patch list
00668     // - this makes it easier to transfer them to an adjacent patch if reqd
00669     {
00670         IFstream is(inputName);
00671 
00672         if (is.good())
00673         {
00674             readHeader(is, fileSignature);
00675 
00676             while ((is >> lineLabel).good())
00677             {
00678                 nFaces++;
00679                 is  >> starCellId
00680                     >> cellFaceId
00681                     >> starRegion
00682                     >> configNumber
00683                     >> patchType;
00684 
00685                 // Build translation table to convert star patch to foam patch
00686                 label patchLabel = mapToFoamPatchId[starRegion];
00687                 if (patchLabel == -1)
00688                 {
00689                     patchLabel = nPatches;
00690                     mapToFoamPatchId[starRegion] = patchLabel;
00691                     origRegion[patchLabel] = starRegion;
00692                     patchTypes_[patchLabel] = patchType;
00693 
00694                     maxId = max(maxId, starRegion);
00695 
00696                     if (patchType == "BAFF")    // should actually be case-insensitive
00697                     {
00698                         nBafflePatches++;
00699                     }
00700                     nPatches++;
00701                 }
00702 
00703                 nPatchFaces[patchLabel]++;
00704             }
00705 
00706             if (nPatches == 0)
00707             {
00708                 Info<< "No boundary faces in file " << inputName << endl;
00709             }
00710         }
00711         else
00712         {
00713             Info<< "Could not read boundary file " << inputName << endl;
00714         }
00715     }
00716 
00717     // keep empty patch region in reserve
00718     nPatches++;
00719     Info<< "Number of patches = " << nPatches
00720         << " (including extra for missing)" << endl;
00721 
00722     // resize
00723     origRegion.setSize(nPatches);
00724     patchTypes_.setSize(nPatches);
00725     patchNames_.setSize(nPatches);
00726     nPatchFaces.setSize(nPatches);
00727 
00728     // add our empty patch
00729     origRegion[nPatches-1] = 0;
00730     nPatchFaces[nPatches-1] = 0;
00731     patchTypes_[nPatches-1] = "none";
00732 
00733     // create names
00734     // - use 'Label' entry from "constant/boundaryRegion" dictionary
00735     forAll(patchTypes_, patchI)
00736     {
00737         bool foundName = false, foundType = false;
00738 
00739         Map<dictionary>::const_iterator
00740             iter = boundaryRegion_.find(origRegion[patchI]);
00741 
00742         if
00743         (
00744             iter != boundaryRegion_.end()
00745         )
00746         {
00747             foundType = iter().readIfPresent
00748             (
00749                 "BoundaryType",
00750                 patchTypes_[patchI]
00751             );
00752 
00753             foundName = iter().readIfPresent
00754             (
00755                 "Label",
00756                 patchNames_[patchI]
00757             );
00758         }
00759 
00760         // consistent names, in long form and in lowercase
00761         if (!foundType)
00762         {
00763             // transform
00764             forAllIter(string, patchTypes_[patchI], i)
00765             {
00766                 *i = tolower(*i);
00767             }
00768 
00769             if (patchTypes_[patchI] == "symp")
00770             {
00771                 patchTypes_[patchI] = "symplane";
00772             }
00773             else if (patchTypes_[patchI] == "cycl")
00774             {
00775                 patchTypes_[patchI] = "cyclic";
00776             }
00777             else if (patchTypes_[patchI] == "baff")
00778             {
00779                 patchTypes_[patchI] = "baffle";
00780             }
00781             else if (patchTypes_[patchI] == "moni")
00782             {
00783                 patchTypes_[patchI] = "monitoring";
00784             }
00785         }
00786 
00787         // create a name if needed
00788         if (!foundName)
00789         {
00790             patchNames_[patchI] =
00791                 patchTypes_[patchI] + "_" + name(origRegion[patchI]);
00792         }
00793     }
00794 
00795     // enforce name "Default_Boundary_Region"
00796     patchNames_[nPatches-1] = defaultBoundaryName;
00797 
00798     // sort according to ascending region numbers, but leave
00799     // Default_Boundary_Region as the final patch
00800     {
00801         labelList sortedIndices;
00802         sortedOrder(SubList<label>(origRegion, nPatches-1), sortedIndices);
00803 
00804         labelList oldToNew = identity(nPatches);
00805         forAll(sortedIndices, i)
00806         {
00807             oldToNew[sortedIndices[i]] = i;
00808         }
00809 
00810         inplaceReorder(oldToNew, origRegion);
00811         inplaceReorder(oldToNew, patchTypes_);
00812         inplaceReorder(oldToNew, patchNames_);
00813         inplaceReorder(oldToNew, nPatchFaces);
00814     }
00815 
00816     // re-sort to have baffles near the end
00817     nBafflePatches = 1;
00818     if (nBafflePatches)
00819     {
00820         labelList oldToNew = identity(nPatches);
00821         label newIndex = 0;
00822         label baffleIndex = (nPatches-1 - nBafflePatches);
00823 
00824         for (label i=0; i < oldToNew.size()-1; ++i)
00825         {
00826             if (patchTypes_[i] == "baffle")
00827             {
00828                 oldToNew[i] = baffleIndex++;
00829             }
00830             else
00831             {
00832                 oldToNew[i] = newIndex++;
00833             }
00834         }
00835 
00836         inplaceReorder(oldToNew, origRegion);
00837         inplaceReorder(oldToNew, patchTypes_);
00838         inplaceReorder(oldToNew, patchNames_);
00839         inplaceReorder(oldToNew, nPatchFaces);
00840     }
00841 
00842     mapToFoamPatchId.setSize(maxId+1, -1);
00843     forAll(origRegion, patchI)
00844     {
00845         mapToFoamPatchId[origRegion[patchI]] = patchI;
00846     }
00847 
00848     boundaryIds_.setSize(nPatches);
00849     forAll(boundaryIds_, patchI)
00850     {
00851         boundaryIds_[patchI].setSize(nPatchFaces[patchI]);
00852         nPatchFaces[patchI] = 0;
00853     }
00854 
00855 
00856     // Pass 2:
00857     //
00858     if (nPatches > 1 && mapToFoamCellId_.size() > 1)
00859     {
00860         IFstream is(inputName);
00861         readHeader(is, fileSignature);
00862 
00863         while ((is >> lineLabel).good())
00864         {
00865             is
00866                 >> starCellId
00867                 >> cellFaceId
00868                 >> starRegion
00869                 >> configNumber
00870                 >> patchType;
00871 
00872             label patchI = mapToFoamPatchId[starRegion];
00873 
00874             // zero-based indexing
00875             cellFaceId--;
00876 
00877             label cellId = -1;
00878 
00879             // convert to FOAM cell number
00880             if (starCellId < mapToFoamCellId_.size())
00881             {
00882                 cellId = mapToFoamCellId_[starCellId];
00883             }
00884 
00885             if (cellId < 0)
00886             {
00887                 Info
00888                     << "Boundaries inconsistent with cell file. "
00889                     << "Star cell " << starCellId << " does not exist"
00890                     << endl;
00891             }
00892             else
00893             {
00894                 // restrict lookup to volume cells (no baffles)
00895                 if (cellId < cellShapes_.size())
00896                 {
00897                     label index = cellShapes_[cellId].model().index();
00898                     if (faceLookupIndex.found(index))
00899                     {
00900                         index = faceLookupIndex[index];
00901                         cellFaceId = starToFoamFaceAddr[index][cellFaceId];
00902                     }
00903                 }
00904                 else
00905                 {
00906                     // we currently use cellId >= nCells to tag baffles,
00907                     // we can also use a negative face number
00908                     cellFaceId = -1;
00909                 }
00910 
00911                 boundaryIds_[patchI][nPatchFaces[patchI]] =
00912                     cellFaceIdentifier(cellId, cellFaceId);
00913 
00914 #ifdef DEBUG_BOUNDARY
00915                 Info<< "bnd " << cellId << " " << cellFaceId << endl;
00916 #endif
00917                 // increment counter of faces in current patch
00918                 nPatchFaces[patchI]++;
00919             }
00920         }
00921     }
00922 
00923     // retain original information in patchPhysicalTypes_ - overwrite latter
00924     patchPhysicalTypes_.setSize(patchTypes_.size());
00925 
00926 
00927     forAll(boundaryIds_, patchI)
00928     {
00929         // resize - avoid invalid boundaries
00930         if (nPatchFaces[patchI] < boundaryIds_[patchI].size())
00931         {
00932             boundaryIds_[patchI].setSize(nPatchFaces[patchI]);
00933         }
00934 
00935         word origType = patchTypes_[patchI];
00936         patchPhysicalTypes_[patchI] = origType;
00937 
00938         if (origType == "symplane")
00939         {
00940             patchTypes_[patchI] = symmetryPolyPatch::typeName;
00941             patchPhysicalTypes_[patchI] = patchTypes_[patchI];
00942         }
00943         else if (origType == "wall")
00944         {
00945             patchTypes_[patchI] = wallPolyPatch::typeName;
00946             patchPhysicalTypes_[patchI] = patchTypes_[patchI];
00947         }
00948         else if (origType == "cyclic")
00949         {
00950             // incorrect. should be cyclicPatch but this
00951             // requires info on connected faces.
00952             patchTypes_[patchI] = cyclicPolyPatch::typeName;
00953             patchPhysicalTypes_[patchI] = patchTypes_[patchI];
00954         }
00955         else if (origType == "baffle")
00956         {
00957             // incorrect. tag the patch until we get proper support.
00958             // set physical type to a canonical "baffle"
00959             patchTypes_[patchI] = emptyPolyPatch::typeName;
00960             patchPhysicalTypes_[patchI] = "baffle";
00961         }
00962         else
00963         {
00964             patchTypes_[patchI] = polyPatch::typeName;
00965         }
00966 
00967         Info<< "patch " << patchI
00968             << " (region " << origRegion[patchI]
00969             << ": " << origType << ") type: '" << patchTypes_[patchI]
00970             << "' name: " << patchNames_[patchI] << endl;
00971     }
00972 
00973     // cleanup
00974     mapToFoamCellId_.clear();
00975     cellShapes_.clear();
00976 }
00977 
00978 
00979 //
00980 // remove unused points
00981 //
00982 void Foam::meshReaders::STARCD::cullPoints()
00983 {
00984     label nPoints = points_.size();
00985     labelList oldToNew(nPoints, -1);
00986 
00987     // loop through cell faces and note which points are being used
00988     forAll(cellFaces_, cellI)
00989     {
00990         const faceList& faces = cellFaces_[cellI];
00991         forAll(faces, i)
00992         {
00993             const labelList& labels = faces[i];
00994             forAll(labels, j)
00995             {
00996                 oldToNew[labels[j]]++;
00997             }
00998         }
00999     }
01000 
01001     // the new ordering and the count of unused points
01002     label pointI = 0;
01003     forAll(oldToNew, i)
01004     {
01005         if (oldToNew[i] >= 0)
01006         {
01007             oldToNew[i] = pointI++;
01008         }
01009     }
01010 
01011     // report unused points
01012     if (nPoints > pointI)
01013     {
01014         Info<< "Unused    points  = " << (nPoints - pointI) << endl;
01015         nPoints = pointI;
01016 
01017         // adjust points and truncate
01018         inplaceReorder(oldToNew, points_);
01019         points_.setSize(nPoints);
01020 
01021         // adjust cellFaces - with mesh shapes this might be faster
01022         forAll(cellFaces_, cellI)
01023         {
01024             faceList& faces = cellFaces_[cellI];
01025             forAll(faces, i)
01026             {
01027                 inplaceRenumber(oldToNew, faces[i]);
01028             }
01029         }
01030 
01031         // adjust baffles
01032         forAll(baffleFaces_, faceI)
01033         {
01034             inplaceRenumber(oldToNew, baffleFaces_[faceI]);
01035         }
01036     }
01037 }
01038 
01039 
01040 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
01041 
01042 bool Foam::meshReaders::STARCD::readGeometry(const scalar scaleFactor)
01043 {
01044     // Info<< "called meshReaders::STARCD::readGeometry" << endl;
01045 
01046     readPoints(geometryFile_ + ".vrt", scaleFactor);
01047     readCells(geometryFile_ + ".cel");
01048     cullPoints();
01049     readBoundary(geometryFile_ + ".bnd");
01050 
01051     return true;
01052 }
01053 
01054 
01055 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
01056 
01057 Foam::meshReaders::STARCD::STARCD
01058 (
01059     const fileName& prefix,
01060     const objectRegistry& registry,
01061     const scalar scaleFactor
01062 )
01063 :
01064     meshReader(prefix, scaleFactor),
01065     cellShapes_(0),
01066     mapToFoamPointId_(0),
01067     mapToFoamCellId_(0)
01068 {
01069     readAux(registry);
01070 }
01071 
01072 
01073 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
01074 
01075 Foam::meshReaders::STARCD::~STARCD()
01076 {}
01077 
01078 
01079 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines