src/surfMesh/surfaceFormats/starcd/STARCDsurfaceFormat.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2009 OpenCFD Ltd.
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "STARCDsurfaceFormat.H"
  28 #include "ListOps.H"
  29
  30 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  31
  32 template<class Face>
  33 inline void Foam::fileFormats::STARCDsurfaceFormat<Face>::writeShell
  34 (
  35     Ostream& os,
  36     const Face& f,
  37     const label cellId,
  38     const label cellTableId
  39 )
  40 {
  41     os  << cellId                    // includes 1 offset
  42         << ' ' << starcdShellShape_  // 3(shell) shape
  43         << ' ' << f.size()
  44         << ' ' << cellTableId
  45         << ' ' << starcdShellType_;  // 4(shell)
  46
  47     // primitives have <= 8 vertices, but prevent overrun anyhow
  48     // indent following lines for ease of reading
  49     label count = 0;
  50     forAll(f, fp)
  51     {
  52         if ((count % 8) == 0)
  53         {
  54             os  << nl << "  " << cellId;
  55         }
  56         os  << ' ' << f[fp] + 1;
  57         count++;
  58     }
  59     os  << endl;
  60 }
  61
  62
  63 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  64
  65 template<class Face>
  66 Foam::fileFormats::STARCDsurfaceFormat<Face>::STARCDsurfaceFormat
  67 (
  68     const fileName& filename
  69 )
  70 {
  71     read(filename);
  72 }
  73
  74
  75 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  76
  77 template<class Face>
  78 bool Foam::fileFormats::STARCDsurfaceFormat<Face>::read
  79 (
  80     const fileName& filename
  81 )
  82 {
  83     const bool mustTriangulate = this->isTri();
  84     this->clear();
  85
  86     fileName baseName = filename.lessExt();
  87
  88     // STAR-CD index of points
  89     List<label> pointId;
  90
  91     // read points from .vrt file
  92     readPoints
  93     (
  94         IFstream(baseName + ".vrt")(),
  95         this->storedPoints(),
  96         pointId
  97     );
  98
  99     // Build inverse mapping (STAR-CD pointId -> index)
 100     Map<label> mapPointId(2*pointId.size());
 101     forAll(pointId, i)
 102     {
 103         mapPointId.insert(pointId[i], i);
 104     }
 105     pointId.clear();
 106
 107     //
 108     // read .cel file
 109     // ~~~~~~~~~~~~~~
 110     IFstream is(baseName + ".cel");
 111     if (!is.good())
 112     {
 113         FatalErrorIn
 114         (
 115             "fileFormats::STARCDsurfaceFormat::read(const fileName&)"
 116         )
 117             << "Cannot read file " << is.name()
 118             << exit(FatalError);
 119     }
 120
 121     readHeader(is, "PROSTAR_CELL");
 122
 123     DynamicList<Face>  dynFaces;
 124     DynamicList<label> dynZones;
 125     DynamicList<word>  dynNames;
 126     DynamicList<label> dynSizes;
 127     Map<label> lookup;
 128
 129     // assume the cellTableIds are not intermixed
 130     bool sorted = true;
 131     label zoneI = 0;
 132
 133     label lineLabel, shapeId, nLabels, cellTableId, typeId;
 134     DynamicList<label> vertexLabels(64);
 135
 136     while ((is >> lineLabel).good())
 137     {
 138         is >> shapeId >> nLabels >> cellTableId >> typeId;
 139
 140         vertexLabels.clear();
 141         vertexLabels.reserve(nLabels);
 142
 143         // read indices - max 8 per line
 144         for (label i = 0; i < nLabels; ++i)
 145         {
 146             label vrtId;
 147             if ((i % 8) == 0)
 148             {
 149                is >> lineLabel;
 150             }
 151             is >> vrtId;
 152
 153             // convert original vertex id to point label
 154             vertexLabels.append(mapPointId[vrtId]);
 155         }
 156
 157         if (typeId == starcdShellType_)
 158         {
 159             // Convert groupID into zoneID
 160             Map<label>::const_iterator fnd = lookup.find(cellTableId);
 161             if (fnd != lookup.end())
 162             {
 163                 if (zoneI != fnd())
 164                 {
 165                     // cellTableIds are intermixed
 166                     sorted = false;
 167                 }
 168                 zoneI = fnd();
 169             }
 170             else
 171             {
 172                 zoneI = dynSizes.size();
 173                 lookup.insert(cellTableId, zoneI);
 174                 dynNames.append(word("cellTable_") + ::Foam::name(zoneI));
 175                 dynSizes.append(0);
 176             }
 177
 178             SubList<label> vertices(vertexLabels, vertexLabels.size());
 179             if (mustTriangulate && nLabels > 3)
 180             {
 181                 face f(vertices);
 182
 183                 faceList triFaces(f.nTriangles());
 184                 label nTri = 0;
 185                 f.triangles(this->points(), nTri, triFaces);
 186
 187                 forAll(triFaces, faceI)
 188                 {
 189                     // a triangular face, but not yet a triFace
 190                     dynFaces.append
 191                     (
 192                         triFace
 193                         (
 194                             static_cast<UList<label>&>(triFaces[faceI])
 195                         )
 196                     );
 197                     dynZones.append(zoneI);
 198                     dynSizes[zoneI]++;
 199                 }
 200             }
 201             else
 202             {
 203                 dynFaces.append(Face(vertices));
 204                 dynZones.append(zoneI);
 205                 dynSizes[zoneI]++;
 206             }
 207         }
 208     }
 209     mapPointId.clear();
 210
 211     sortFacesAndStore(dynFaces.xfer(), dynZones.xfer(), sorted);
 212
 213     // add zones, culling empty ones
 214     this->addZones(dynSizes, dynNames, true);
 215     return true;
 216 }
 217
 218
 219 template<class Face>
 220 void Foam::fileFormats::STARCDsurfaceFormat<Face>::write
 221 (
 222     const fileName& filename,
 223     const MeshedSurfaceProxy<Face>& surf
 224 )
 225 {
 226     const pointField& pointLst = surf.points();
 227     const List<Face>&  faceLst = surf.faces();
 228     const List<label>& faceMap = surf.faceMap();
 229
 230     const List<surfZone>& zones =
 231     (
 232         surf.surfZones().size() > 1
 233       ? surf.surfZones()
 234       : oneZone(faceLst)
 235     );
 236
 237     const bool useFaceMap = (surf.useFaceMap() && zones.size() > 1);
 238
 239
 240     fileName baseName = filename.lessExt();
 241
 242     writePoints(OFstream(baseName + ".vrt")(), pointLst);
 243     OFstream os(baseName + ".cel");
 244     writeHeader(os, "CELL");
 245
 246     label faceIndex = 0;
 247     forAll(zones, zoneI)
 248     {
 249         const surfZone& zone = zones[zoneI];
 250
 251         if (useFaceMap)
 252         {
 253             forAll(zone, localFaceI)
 254             {
 255                 const Face& f = faceLst[faceMap[faceIndex++]];
 256                 writeShell(os, f, faceIndex, zoneI + 1);
 257             }
 258         }
 259         else
 260         {
 261             forAll(zone, localFaceI)
 262             {
 263                 const Face& f = faceLst[faceIndex++];
 264                 writeShell(os, f, faceIndex, zoneI + 1);
 265             }
 266         }
 267     }
 268
 269     // write simple .inp file
 270     writeCase
 271     (
 272         OFstream(baseName + ".inp")(),
 273         pointLst,
 274         faceLst.size(),
 275         zones
 276     );
 277 }
 278
 279
 280 // ************************************************************************* //