Line data Source code
1 : /*
2 : * Copyright (c) 2014-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Sebastian Reiter
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : #ifndef __H__UG_file_io_vtu_impl
34 : #define __H__UG_file_io_vtu_impl
35 :
36 : #include <sstream>
37 : #include <cstring>
38 : #include "lib_grid/algorithms/debug_util.h"
39 : #include "lib_grid/callbacks/subset_callbacks.h"
40 :
41 : namespace ug{
42 :
43 :
44 : template <class TAPosition>
45 0 : bool LoadGridFromVTU(Grid& grid, ISubsetHandler& sh, const char* filename,
46 : TAPosition& aPos)
47 : {
48 0 : GridReaderVTU vtuReader;
49 0 : if(!vtuReader.parse_file(filename)){
50 0 : UG_LOG("ERROR in LoadGridFromVTU: File not found: " << filename << std::endl);
51 : return false;
52 : }
53 :
54 0 : if(vtuReader.num_grids() < 1){
55 : UG_LOG("ERROR in LoadGridFromVTU: File contains no grid.\n");
56 : return false;
57 : }
58 :
59 0 : vtuReader.grid(grid, 0, aPos);
60 :
61 0 : vtuReader.subset_handler(sh, 0, 0);
62 :
63 : return true;
64 0 : }
65 :
66 : template <class TAPosition>
67 0 : bool SaveGridToVTU(Grid& grid, ISubsetHandler* psh, const char* filename,
68 : TAPosition& aPos)
69 : {
70 0 : GridWriterVTU vtuWriter;
71 0 : std::ofstream out(filename);
72 :
73 0 : vtuWriter.set_stream(&out);
74 :
75 0 : vtuWriter.new_piece(grid, psh, aPos);
76 :
77 0 : vtuWriter.finish();
78 :
79 0 : return true;
80 0 : }
81 :
82 :
83 : ////////////////////////////////////////////////////////////////////////////////
84 : ////////////////////////////////////////////////////////////////////////////////
85 : // GridWriterVTU
86 : ////////////////////////////////////////////////////////////////////////////////
87 :
88 0 : std::ostream& GridWriterVTU::
89 : out_stream()
90 : {
91 0 : UG_COND_THROW(!m_pout, "Invalid ostream specified fo GridWriterVTU!");
92 0 : return *m_pout;
93 : }
94 :
95 : template <class TPositionAttachment>
96 0 : bool GridWriterVTU::
97 : new_piece(Grid& grid, ISubsetHandler* psh, TPositionAttachment& aPos)
98 : {
99 : using namespace std;
100 :
101 0 : if(m_pieceMode == OPEN)
102 0 : end_piece();
103 :
104 0 : m_pieceMode = OPEN;
105 0 : m_pointDataMode = NONE;
106 0 : m_cellDataMode = NONE;
107 :
108 0 : ostream& out = out_stream();
109 :
110 0 : m_curGrid = &grid;
111 0 : m_curSH = psh;
112 :
113 : m_pieceSubsetHandlers.clear();
114 : m_cells.clear();
115 :
116 : // if a subset-handler is specified, all grid elements which are assigned to
117 : // subsets are considered as cells.
118 : // If not, only elements of highest dimension are considered as cells.
119 :
120 0 : if(psh){
121 0 : add_subset_handler(*psh, string("regions"));
122 0 : collect_cells<Vertex>(m_cells, grid, IsNotInSubset(*psh, -1));
123 0 : collect_cells<Edge>(m_cells, grid, IsNotInSubset(*psh, -1));
124 0 : collect_cells<Face>(m_cells, grid, IsNotInSubset(*psh, -1));
125 0 : collect_cells<Volume>(m_cells, grid, IsNotInSubset(*psh, -1));
126 : }
127 : else{
128 0 : if(grid.num<Volume>() > 0)
129 0 : collect_cells<Volume>(m_cells, grid, ConsiderAll());
130 0 : else if(grid.num<Face>() > 0)
131 0 : collect_cells<Face>(m_cells, grid, ConsiderAll());
132 0 : else if(grid.num<Edge>() > 0)
133 0 : collect_cells<Edge>(m_cells, grid, ConsiderAll());
134 : else
135 0 : collect_cells<Vertex>(m_cells, grid, ConsiderAll());
136 : }
137 :
138 : out << " <Piece NumberOfPoints=\"" << grid.num<Vertex>() << "\""
139 : << " NumberOfCells=\"" << m_cells.size() << "\">" << endl;
140 :
141 :
142 : // write points
143 : out << " <Points>" << endl;
144 0 : write_vector_data<Vertex>(grid, aPos, "");
145 : out << " </Points>" << endl;
146 :
147 : // write cells
148 : // first we'll assign indices to the vertices, which can then be used as
149 : // references by the cells.
150 : AInt aInd;
151 : grid.attach_to_vertices(aInd);
152 : Grid::VertexAttachmentAccessor<AInt> aaInd(grid, aInd);
153 : AssignIndices(grid.begin<Vertex>(), grid.end<Vertex>(), aaInd, 0);
154 :
155 0 : write_cells(m_cells, grid, aaInd);
156 :
157 : grid.detach_from_vertices(aInd);
158 :
159 0 : return true;
160 : }
161 :
162 :
163 : template <class TElem, class TAttachment>
164 0 : void GridWriterVTU::
165 : write_vector_data(Grid& grid,
166 : TAttachment aData,
167 : const char* name,
168 : typename Grid::traits<TElem>::callback consider_elem)
169 : {
170 : using namespace std;
171 0 : ostream& out = out_stream();
172 :
173 : typedef typename TAttachment::ValueType vector_t;
174 : typedef typename Grid::traits<TElem>::iterator iter_t;
175 :
176 : Grid::AttachmentAccessor<TElem, TAttachment> aaData(grid, aData);
177 :
178 0 : write_data_array_header("Float32", "", vector_t::Size);
179 0 : out << " ";
180 :
181 0 : for(iter_t i = grid.begin<TElem>(); i != grid.end<TElem>(); ++i){
182 : vector_t& v = aaData[*i];
183 0 : for(size_t c = 0; c < vector_t::Size; ++c){
184 0 : out << " " << v[c];
185 : }
186 : }
187 : out << endl;
188 0 : write_data_array_footer();
189 0 : }
190 :
191 : template <class TElem>
192 0 : void GridWriterVTU::
193 : collect_cells(std::vector<GridObject*>& cellsOut, Grid& grid,
194 : typename Grid::traits<TElem>::callback consider_elem)
195 : {
196 : typedef typename Grid::traits<TElem>::iterator iter_t;
197 0 : for(iter_t i = grid.begin<TElem>(); i != grid.end<TElem>(); ++i){
198 0 : if(consider_elem(*i))
199 0 : cellsOut.push_back(*i);
200 : }
201 0 : }
202 :
203 :
204 :
205 :
206 : ////////////////////////////////////////////////////////////////////////////////
207 : ////////////////////////////////////////////////////////////////////////////////
208 : // GridReaderVTU
209 : ////////////////////////////////////////////////////////////////////////////////
210 : template <class TPositionAttachment>
211 0 : bool GridReaderVTU::
212 : grid(Grid& gridOut, size_t index, TPositionAttachment& aPos)
213 : {
214 : using namespace rapidxml;
215 : using namespace std;
216 :
217 : // make sure that a node at the given index exists
218 0 : if(num_grids() <= index){
219 : UG_LOG(" GridReaderVTU::read: bad grid index!\n");
220 0 : return false;
221 : }
222 :
223 : Grid& grid = gridOut;
224 :
225 : // Since we have to create all elements in the correct order and
226 : // since we have to make sure that no elements are created in between,
227 : // we'll first disable all grid-options and reenable them later on
228 0 : uint gridopts = grid.get_options();
229 0 : grid.set_options(GRIDOPT_NONE);
230 :
231 : // access node data
232 0 : if(!grid.has_vertex_attachment(aPos)){
233 0 : grid.attach_to_vertices(aPos);
234 : }
235 :
236 : Grid::VertexAttachmentAccessor<TPositionAttachment> aaPos(grid, aPos);
237 :
238 : // store the grid in the grid-vector and assign indices to the vertices
239 0 : m_entries[index].grid = &grid;
240 :
241 : // get the grid-node and the vertex-vector
242 0 : xml_node<>* gridNode = m_entries[index].node;
243 0 : vector<Vertex*>& vertices = m_entries[index].vertices;
244 0 : vector<GridObject*>& cells = m_entries[index].cells;
245 :
246 : // iterate over all pieces
247 : xml_node<>* pieceNode = gridNode;
248 : // first we'll create all points and cells, then we'll parse point- and cell-data
249 0 : xml_node<>* pointsNode = pieceNode->first_node("Points");
250 0 : UG_COND_THROW(pointsNode == NULL, "Missing Points node in UnstructuredGrid node!")
251 :
252 : size_t vrtOffset = vertices.size();
253 0 : create_vertices(vertices, grid, pointsNode, aaPos);
254 :
255 0 : xml_node<>* cellsNode = pieceNode->first_node("Cells");
256 0 : UG_COND_THROW(cellsNode == NULL, "Missing Cells node in UnstructuredGrid node!")
257 0 : create_cells(cells, grid, cellsNode, vertices, vrtOffset);
258 :
259 : // reenable the grids options.
260 0 : grid.set_options(gridopts);
261 :
262 : return true;
263 : }
264 :
265 : template <class TAAPos>
266 0 : bool GridReaderVTU::
267 : create_vertices(std::vector<Vertex*>& vrtsOut, Grid& grid,
268 : rapidxml::xml_node<>* vrtNode, TAAPos aaPos)
269 : {
270 : using namespace rapidxml;
271 : using namespace std;
272 :
273 0 : xml_node<>* dataNode = vrtNode->first_node("DataArray");
274 0 : UG_COND_THROW(!dataNode, "Missing DataArray node in Points node");
275 :
276 0 : int numSrcCoords = -1;
277 0 : xml_attribute<>* attrib = dataNode->first_attribute("NumberOfComponents");
278 0 : if(attrib)
279 0 : numSrcCoords = atoi(attrib->value());
280 :
281 0 : int numDestCoords = (int)TAAPos::ValueType::Size;
282 :
283 : assert(numDestCoords > 0 && "bad position attachment type");
284 :
285 0 : if(numSrcCoords < 1 || numDestCoords < 1)
286 : return false;
287 :
288 : // create a buffer with which we can access the data
289 0 : string str(dataNode->value(), dataNode->value_size());
290 0 : stringstream ss(str, ios_base::in);
291 :
292 : // if numDestCoords == numSrcCoords parsing will be faster
293 0 : if(numSrcCoords == numDestCoords){
294 0 : while(!ss.eof()){
295 : // read the data
296 : typename TAAPos::ValueType v;
297 :
298 0 : for(int i = 0; i < numSrcCoords; ++i)
299 0 : ss >> v[i];
300 :
301 : // make sure that everything went right
302 0 : if(ss.fail())
303 : break;
304 :
305 : // create a new vertex
306 0 : RegularVertex* vrt = *grid.create<RegularVertex>();
307 0 : vrtsOut.push_back(vrt);
308 :
309 : // set the coordinates
310 : aaPos[vrt] = v;
311 : }
312 : }
313 : else{
314 : // we have to be careful with reading.
315 : // if numDestCoords < numSrcCoords we'll ignore some coords,
316 : // in the other case we'll add some 0's.
317 0 : int minNumCoords = min(numSrcCoords, numDestCoords);
318 0 : typename TAAPos::ValueType::value_type dummy = 0;
319 :
320 0 : while(!ss.eof()){
321 : // read the data
322 : typename TAAPos::ValueType v;
323 :
324 : int iMin;
325 0 : for(iMin = 0; iMin < minNumCoords; ++iMin)
326 0 : ss >> v[iMin];
327 :
328 : // ignore unused entries in the input buffer
329 0 : for(int i = iMin; i < numSrcCoords; ++i)
330 : ss >> dummy;
331 :
332 : // add 0's to the vector
333 0 : for(int i = iMin; i < numDestCoords; ++i)
334 0 : v[i] = 0;
335 :
336 : // make sure that everything went right
337 0 : if(ss.fail()){
338 : UG_LOG("GridReaderVTU::create_vertices: Failed to read vertex.\n");
339 0 : return false;
340 : }
341 :
342 : // create a new vertex
343 0 : RegularVertex* vrt = *grid.create<RegularVertex>();
344 0 : vrtsOut.push_back(vrt);
345 :
346 : // set the coordinates
347 : aaPos[vrt] = v;
348 : }
349 : }
350 :
351 : return true;
352 0 : }
353 :
354 :
355 : template <class T>
356 0 : void GridReaderVTU::
357 : read_scalar_data(std::vector<T>& dataOut,
358 : rapidxml::xml_node<>* dataNode,
359 : bool clearData)
360 : {
361 : using namespace std;
362 :
363 0 : if(clearData)
364 : dataOut.clear();
365 :
366 : // create a buffer with which we can access the data
367 0 : string str(dataNode->value(), dataNode->value_size());
368 0 : stringstream ss(str, ios_base::in);
369 :
370 0 : while(!ss.eof()){
371 : // read the data
372 : T d;
373 0 : ss >> d;
374 :
375 : // make sure that everything went right
376 0 : if(ss.fail())
377 : break;
378 :
379 0 : dataOut.push_back(d);
380 : }
381 0 : }
382 :
383 : template <class T>
384 0 : void GridReaderVTU::
385 : check_indices(std::vector<T>& inds, size_t first, size_t num, size_t validSize)
386 : {
387 0 : UG_COND_THROW(first + num > inds.size(),
388 : "Bad index range encountered during parsing of " << m_filename);
389 :
390 0 : for(size_t i = first; i < num; ++i){
391 0 : UG_COND_THROW(inds[i] >= validSize,
392 : "Bad index encountered during parsing of " << m_filename);
393 : }
394 0 : }
395 :
396 : }// end of namespace
397 :
398 : #endif //__H__file_io_vtu_impl
|
