18 #ifndef DUNE_GRIDGLUE_ADAPTER_GRIDGLUEVTKWRITER_HH 19 #define DUNE_GRIDGLUE_ADAPTER_GRIDGLUEVTKWRITER_HH 26 #include <dune/common/classname.hh> 27 #include <dune/common/typetraits.hh> 28 #include <dune/geometry/type.hh> 29 #include <dune/geometry/referenceelements.hh> 42 template <
class Glue,
int s
ide>
43 static void writeExtractedPart(
const Glue& glue,
const std::string& filename)
45 static_assert((side==0 || side==1),
"'side' can only be 0 or 1");
49 fgrid.open(filename.c_str());
51 typedef typename Dune::conditional<(side==0), typename Glue::Grid0View, typename Glue::Grid1View>::type GridView;
52 typedef typename Dune::conditional<(side==0), typename Glue::Grid0Patch, typename Glue::Grid1Patch>::type
Extractor;
54 typedef typename GridView::ctype ctype;
56 const int domdimw = GridView::dimensionworld;
60 std::string coordinatePadding;
61 for (
int i=domdimw; i<3; i++)
62 coordinatePadding +=
" 0";
64 fgrid <<
"# vtk DataFile Version 2.0\nFilename: " << filename <<
"\nASCII" << std::endl;
68 std::vector<typename Extractor::Coords> coords;
69 glue.template patch<side>().getCoords(coords);
71 fgrid << ((patchDim==3) ?
"DATASET UNSTRUCTURED_GRID" :
"DATASET POLYDATA") << std::endl;
72 fgrid <<
"POINTS " << coords.size() <<
" " << Dune::className<ctype>() << std::endl;
74 for (
size_t i=0; i<coords.size(); i++)
75 fgrid << coords[i] << coordinatePadding << std::endl;
82 std::vector<typename Extractor::VertexVector> faces;
83 std::vector<Dune::GeometryType> geometryTypes;
84 glue.template patch<side>().getFaces(faces);
85 glue.template patch<side>().getGeometryTypes(geometryTypes);
87 unsigned int faceCornerCount = 0;
88 for (
size_t i=0; i<faces.size(); i++)
89 faceCornerCount += faces[i].size();
91 fgrid << ((patchDim==3) ?
"CELLS " :
"POLYGONS ")
92 << geometryTypes.size() <<
" " << geometryTypes.size() + faceCornerCount << std::endl;
94 for (
size_t i=0; i<faces.size(); i++) {
96 fgrid << faces[i].size();
100 if (geometryTypes[i].isSimplex()) {
101 for (
int j=0; j<patchDim+1; j++)
102 fgrid <<
" " << faces[i][j];
104 }
else if (geometryTypes[i].isQuadrilateral()) {
105 fgrid <<
" " << faces[i][0] <<
" " << faces[i][1]
106 <<
" " << faces[i][3] <<
" " << faces[i][2];
108 }
else if (geometryTypes[i].isPyramid()) {
109 fgrid <<
" " << faces[i][0] <<
" " << faces[i][1]
110 <<
" " << faces[i][3] <<
" " << faces[i][2] <<
" " << faces[i][4];
112 }
else if (geometryTypes[i].isPrism()) {
113 fgrid <<
" " << faces[i][0] <<
" " << faces[i][2] <<
" " << faces[i][1]
114 <<
" " << faces[i][3] <<
" " << faces[i][5] <<
" " << faces[i][4];
116 }
else if (geometryTypes[i].isHexahedron()) {
117 fgrid <<
" " << faces[i][0] <<
" " << faces[i][1]
118 <<
" " << faces[i][3] <<
" " << faces[i][2]
119 <<
" " << faces[i][4] <<
" " << faces[i][5]
120 <<
" " << faces[i][7] <<
" " << faces[i][6];
123 DUNE_THROW(Dune::NotImplemented,
"Geometry type " << geometryTypes[i] <<
" not supported yet");
134 fgrid <<
"CELL_TYPES " << geometryTypes.size() << std::endl;
136 for (
size_t i=0; i<geometryTypes.size(); i++) {
137 if (geometryTypes[i].isSimplex())
138 fgrid <<
"10" << std::endl;
139 else if (geometryTypes[i].isHexahedron())
140 fgrid <<
"12" << std::endl;
141 else if (geometryTypes[i].isPrism())
142 fgrid <<
"13" << std::endl;
143 else if (geometryTypes[i].isPyramid())
144 fgrid <<
"14" << std::endl;
146 DUNE_THROW(Dune::NotImplemented,
"Geometry type " << geometryTypes[i] <<
" not supported yet");
155 ctype accum = 0.0, delta = 1.0 / (ctype) (gridSubEntityData.size()-1);
157 fgrid <<
"CELL_DATA " << gridSubEntityData.size() << std::endl;
158 fgrid <<
"SCALARS property_coding " << Dune::className<ctype>() <<
" 1" << std::endl;
159 fgrid <<
"LOOKUP_TABLE default" << std::endl;
161 for (
typename GridSubEntityData::const_iterator sEIt = gridSubEntityData.begin();
162 sEIt != gridSubEntityData.end();
163 ++sEIt, accum += delta)
166 fgrid << accum << std::endl;
176 template <
class Glue,
int s
ide>
177 static void writeIntersections(
const Glue& glue,
const std::string& filename)
179 static_assert((side==0 || side==1),
"'side' can only be 0 or 1");
181 std::ofstream fmerged;
183 fmerged.open(filename.c_str());
185 typedef typename Dune::conditional<(side==0), typename Glue::Grid0View, typename Glue::Grid1View>::type GridView;
186 typedef typename Dune::conditional<(side==0),
187 typename Glue::Grid0IntersectionIterator,
188 typename Glue::Grid1IntersectionIterator>::type RemoteIntersectionIterator;
190 typedef typename GridView::ctype ctype;
192 const int domdimw = GridView::dimensionworld;
193 const int intersectionDim = Glue::Intersection::mydim;
196 std::string coordinatePadding;
197 for (
int i=domdimw; i<3; i++)
198 coordinatePadding +=
" 0";
200 int overlaps = glue.size();
204 typedef typename Glue::Grid0Patch
Extractor;
205 std::vector<typename Extractor::Coords> coords;
206 glue.template patch<side>().getCoords(coords);
209 fmerged <<
"# vtk DataFile Version 2.0\nFilename: " << filename <<
"\nASCII" << std::endl;
210 fmerged << ((intersectionDim==3) ?
"DATASET UNSTRUCTURED_GRID" :
"DATASET POLYDATA") << std::endl;
211 fmerged <<
"POINTS " << overlaps*(intersectionDim+1) <<
" " << Dune::className<ctype>() << std::endl;
213 for (RemoteIntersectionIterator isIt = glue.template ibegin<side>();
214 isIt != glue.template iend<side>();
217 for (
int i = 0; i < isIt->geometry().corners(); ++i)
218 fmerged << isIt->geometry().corner(i) << coordinatePadding << std::endl;
224 std::vector<typename Extractor::VertexVector> faces;
225 std::vector<Dune::GeometryType> geometryTypes;
226 glue.template patch<side>().getFaces(faces);
227 glue.template patch<side>().getGeometryTypes(geometryTypes);
229 unsigned int faceCornerCount = 0;
230 for (
size_t i=0; i<faces.size(); i++)
231 faceCornerCount += faces[i].size();
233 int grid0SimplexCorners = intersectionDim+1;
234 fmerged << ((intersectionDim==3) ?
"CELLS " :
"POLYGONS ")
235 << overlaps <<
" " << (grid0SimplexCorners+1)*overlaps << std::endl;
237 for (
int i = 0; i < overlaps; ++i) {
238 fmerged << grid0SimplexCorners;
239 for (
int j=0; j<grid0SimplexCorners; j++)
240 fmerged <<
" " << grid0SimplexCorners*i+j;
241 fmerged << std::endl;
245 if (intersectionDim==3) {
247 fmerged <<
"CELL_TYPES " << overlaps << std::endl;
249 for (
int i = 0; i < overlaps; i++)
250 fmerged <<
"10" << std::endl;
257 ctype accum = 0.0, delta = 1.0 / (ctype) (gridSubEntityData.size()-1);
259 fmerged <<
"CELL_DATA " << overlaps << std::endl;
260 fmerged <<
"SCALARS property_coding " << Dune::className<ctype>() <<
" 1" << std::endl;
261 fmerged <<
"LOOKUP_TABLE default" << std::endl;
263 for (
typename GridSubEntityData::const_iterator sEIt = gridSubEntityData.begin();
264 sEIt != gridSubEntityData.end();
265 ++sEIt, accum += delta)
268 for (
int j = 0; j < sEIt->first.second; ++j)
269 fmerged << accum << std::endl;
276 template<
typename Glue>
277 static void write(
const Glue& glue,
const std::string& filenameTrunk)
281 writeExtractedPart<Glue,0>(glue,
282 filenameTrunk +
"-grid0.vtk");
284 writeIntersections<Glue,0>(glue,
285 filenameTrunk +
"-intersections-grid0.vtk");
288 writeExtractedPart<Glue,1>(glue,
289 filenameTrunk +
"-grid1.vtk");
291 writeIntersections<Glue,1>(glue,
292 filenameTrunk +
"-intersections-grid1.vtk");
301 #endif // DUNE_GRIDGLUE_ADAPTER_GRIDGLUEVTKWRITER_HH static void write(const Glue &glue, const std::string &filenameTrunk)
Definition: gridgluevtkwriter.hh:277
Write remote intersections to a vtk file for debugging purposes.
Definition: gridgluevtkwriter.hh:36
Definition: gridglue.hh:33
Definition: extractor.hh:54
Definition: extractor.hh:53
Provides codimension-independent methods for grid extraction.
Definition: extractor.hh:47