Dune Core Modules (2.8.0)

starcdreader.hh
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_STARCD_READER_HH
4#define DUNE_STARCD_READER_HH
5
7
10#include <iostream>
11#include <fstream>
12#include <memory>
13
14namespace Dune {
15
49 template <class GridType>
51
52 public:
53
63 static std::unique_ptr<GridType> read(const std::string& fileName, bool verbose = true)
64 {
65 // extract the grid dimension
66 const int dim = GridType::dimension;
67
68 // currently only dim = 3 is implemented
69 if (dim != 3)
71 "Reading Star-CD format is not implemented for dimension " << dim);
72
73 // set up the grid factory
75
76 // set the name of the vertex file
77 std::string vertexFileName = fileName + ".vrt";
78
79 // set the vertex input stream
80 std::ifstream vertexFile(vertexFileName.c_str());
81 if (!vertexFile)
82 DUNE_THROW(Dune::IOError, "Could not open " << vertexFileName);
83
84 // read the vertices
85 int dummyIdx;
86 int numberOfVertices = 0;
87 while (vertexFile >> dummyIdx) {
88 numberOfVertices++;
89
91
92 for (int k = 0; k < dim; k++)
93 vertexFile >> position[k];
94
95 factory.insertVertex(position);
96 }
97 if (verbose)
98 std::cout << numberOfVertices << " vertices read." << std::endl;
99
100 // set the name of the element file
101 std::string elementFileName = fileName + ".cel";
102
103 // set the element input stream
104 std::ifstream elementFile(elementFileName.c_str());
105 if (!elementFile)
106 DUNE_THROW(Dune::IOError, "Could not open " << elementFileName);
107
108 // read the elements
109 int numberOfElements = 0;
110 int numberOfSimplices = 0;
111 int numberOfPyramids = 0;
112 int numberOfPrisms = 0;
113 int numberOfCubes = 0;;
114 int maxNumberOfVertices = (int)pow(2, dim);
115 int isVolume = 1;
116 while (elementFile >> dummyIdx) {
117 std::vector<unsigned int> vertices(maxNumberOfVertices);
118 for (int k = 0; k < maxNumberOfVertices; k++)
119 elementFile >> vertices[k];
120
121 int boundaryId;
122 elementFile >> boundaryId;
123
124 int volumeOrSurface[2];
125 elementFile >> volumeOrSurface[0] >> volumeOrSurface[1];
126
127 if (volumeOrSurface[0] == isVolume) {
128 numberOfElements++;
129
130 if (vertices[2] == vertices[3]) { // simplex or prism
131 if (vertices[4] == vertices[5]) { // simplex
132 numberOfSimplices++;
133 std::vector<unsigned int> simplexVertices(4);
134 for (int k = 0; k < 3; k++)
135 simplexVertices[k] = vertices[k] - 1;
136 simplexVertices[3] = vertices[4] - 1;
137 factory.insertElement(Dune::GeometryTypes::tetrahedron, simplexVertices);
138 }
139 else { // prism
140 numberOfPrisms++;
141 std::vector<unsigned int> prismVertices(6);
142 for (int k = 0; k < 3; k++)
143 prismVertices[k] = vertices[k] - 1;
144 for (int k = 3; k < 6; k++)
145 prismVertices[k] = vertices[k+1] - 1;
146 factory.insertElement(Dune::GeometryTypes::prism, prismVertices);
147 }
148 }
149 else { // cube or pyramid
150 if (vertices[4] == vertices[5]) { // pyramid
151 numberOfPyramids++;
152 std::vector<unsigned int> pyramidVertices(5);
153 for (int k = 0; k < 5; k++)
154 pyramidVertices[k] = vertices[k] - 1;
155 factory.insertElement(Dune::GeometryTypes::pyramid, pyramidVertices);
156 }
157 else { // cube
158 numberOfCubes++;
159 std::vector<unsigned int> cubeVertices(8);
160 for (int k = 0; k < 8; k++)
161 cubeVertices[k] = vertices[k] - 1;
162 std::swap(cubeVertices[2], cubeVertices[3]);
163 std::swap(cubeVertices[6], cubeVertices[7]);
164 factory.insertElement(Dune::GeometryTypes::hexahedron, cubeVertices);
165 }
166 }
167 }
168 }
169 if (verbose)
170 std::cout << numberOfElements << " elements read: "
171 << numberOfSimplices << " simplices, " << numberOfPyramids << " pyramids, "
172 << numberOfPrisms << " prisms, " << numberOfCubes << " cubes." << std::endl;
173
174 // finish off the construction of the grid object
175 if (verbose)
176 std::cout << "Starting createGrid() ... " << std::flush;
177
178 return factory.createGrid();
179
180 }
181
182 };
183
184}
185
186#endif
vector space out of a tensor product of fields.
Definition: fvector.hh:95
Provide a generic factory class for unstructured grids.
Definition: gridfactory.hh:312
virtual void insertElement(const GeometryType &type, const std::vector< unsigned int > &vertices)
Insert an element into the coarse grid.
Definition: gridfactory.hh:344
virtual void insertVertex(const FieldVector< ctype, dimworld > &pos)
Insert a vertex into the coarse grid.
Definition: gridfactory.hh:333
virtual std::unique_ptr< GridType > createGrid()
Finalize grid creation and hand over the grid.
Definition: gridfactory.hh:370
Default exception class for I/O errors.
Definition: exceptions.hh:229
Default exception for dummy implementations.
Definition: exceptions.hh:261
File reader for the Star-CD format.
Definition: starcdreader.hh:50
static std::unique_ptr< GridType > read(const std::string &fileName, bool verbose=true)
Read grid from a Star-CD file.
Definition: starcdreader.hh:63
Provide a generic factory class for unstructured grids.
A few common exception classes.
#define DUNE_THROW(E, m)
Definition: exceptions.hh:216
constexpr GeometryType prism
GeometryType representing a 3D prism.
Definition: type.hh:540
constexpr GeometryType hexahedron
GeometryType representing a hexahedron.
Definition: type.hh:546
constexpr GeometryType pyramid
GeometryType representing a 3D pyramid.
Definition: type.hh:534
constexpr GeometryType tetrahedron
GeometryType representing a tetrahedron.
Definition: type.hh:528
Dune namespace.
Definition: alignedallocator.hh:11
A unique label for each type of element that can occur in a grid.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 21, 23:30, 2024)