VtkWriter.cc 12.9 KB
Newer Older
1
#include <stdio.h>
2
3
#include <string>
#include <fstream>
Peter Gottschling's avatar
Peter Gottschling committed
4
#include <cmath>
5
6
7
8
9
10
11
12
13

#include "VtkWriter.h"
#include "DataCollector.h"
#include "DOFVector.h"
#include "SurfaceRegion_ED.h"
#include "ElementRegion_ED.h"

namespace AMDiS { 

14
15
16
17
18
  int VtkWriter::writeFile(const char *name)
  {
    int nVertices = (*dc_)[0]->getNumberVertices();
    int nElements = (*dc_)[0]->getNumberElements();
    int vertices = (*dc_)[0]->getMesh()->getGeo(VERTEX);
19

20
    if ((dim_ == 2) && (degree_ == 2)) {
21
22
      nVertices += (*dc_)[0]->getNumberInterpPoints();
      nElements *= 4;
23
24
25
    } else if ((dim_ == 2) && (degree_ == 3)) {
      nVertices += (*dc_)[0]->getNumberInterpPoints();
      nElements *= 9;
26
27
28
    } else if ((dim_ == 2) && (degree_ == 4)) {
      nVertices += (*dc_)[0]->getNumberInterpPoints();
      nElements *= 16;
29
    }
30
    
Thomas Witkowski's avatar
Thomas Witkowski committed
31
    std::ofstream file;
32
    file.open(name);
33
    
Thomas Witkowski's avatar
Thomas Witkowski committed
34
35
36
    file << "<?xml version=\"1.0\"?>\n";
    file << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n";
    file << "  <UnstructuredGrid>\n";
37
    file << "    <Piece NumberOfPoints=\"" << nVertices 
Thomas Witkowski's avatar
Thomas Witkowski committed
38
39
40
	 << "\" NumberOfCells=\"" << nElements << "\">\n";
    file << "      <Points>\n";
    file << "        <DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\"ascii\">\n";
41
    
42
    writeVertexCoords(file);
43
    
Thomas Witkowski's avatar
Thomas Witkowski committed
44
45
46
47
    file << "        </DataArray>\n";
    file << "      </Points>\n";
    file << "      <Cells>\n";
    file << "        <DataArray type=\"Int32\" Name=\"offsets\">\n";
48
    
49
    for (int i = 0; i < nElements; i++) {
Thomas Witkowski's avatar
Thomas Witkowski committed
50
      file << " " << (i + 1) * vertices << "\n";
51
52
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
53
54
    file << "        </DataArray>\n";
    file << "        <DataArray type=\"UInt8\" Name=\"types\">\n";
55
56

    for (int i = 0; i < nElements; i++) {
57
58
      switch (vertices) {
      case 2:
Thomas Witkowski's avatar
Thomas Witkowski committed
59
	file << " 3\n";
60
61
	break;
      case 3: 
Thomas Witkowski's avatar
Thomas Witkowski committed
62
	file << " 5\n";
63
64
	break;
      case 4:
Thomas Witkowski's avatar
Thomas Witkowski committed
65
	file << " 10\n";
66
67
68
69
70
71
	break;		    
      default:
	break;
      }	   
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
72
73
    file << "        </DataArray>\n";
    file << "        <DataArray type=\"Int32\" Name=\"connectivity\">\n";
74
    
75
76
    writeConnectivity(file);
    
Thomas Witkowski's avatar
Thomas Witkowski committed
77
78
79
    file << "        </DataArray>\n";    
    file << "      </Cells>\n";
    file << "      <PointData>\n";
80
81
    
    for (int i = 0; i < static_cast<int>(dc_->size()); i++) {
82
      file << "        <DataArray type=\"Float32\" Name=\"value" << i 
Thomas Witkowski's avatar
Thomas Witkowski committed
83
	   << "\" format=\"ascii\">\n";
84
      
85
86
      writeVertexValues(file, i);
      
Thomas Witkowski's avatar
Thomas Witkowski committed
87
      file << "        </DataArray>\n";
88
    }
89
    
Thomas Witkowski's avatar
Thomas Witkowski committed
90
91
92
93
    file << "      </PointData>\n";
    file << "    </Piece>\n";
    file << "  </UnstructuredGrid>\n";
    file << "</VTKFile>\n";
94
    
95
    
96
    file.close();
97
    
98
99
    return 0;
  }
100

101

Thomas Witkowski's avatar
Thomas Witkowski committed
102
  void VtkWriter::writeVertexCoords(std::ofstream &file)
103
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
104
105
    DOFVector< std::list<VertexInfo> > *vertexInfos = (*dc_)[0]->getVertexInfos();
    DOFVector< std::list<VertexInfo> >::Iterator it(vertexInfos, USED_DOFS);
106
107
    int counter = 0;

108
    // For all DOFs of vertices, write the coordinates.
109
110
    for (it.reset(); !it.end(); ++it) {
      // for all vertex infos of this DOF
Thomas Witkowski's avatar
Thomas Witkowski committed
111
      std::list<VertexInfo>::iterator it2;
112
113
114
115
116
117
      for (it2 = it->begin(); it2 != it->end(); ++it2) {
	it2->outputIndex = counter++;
	writeCoord(file, it2->coords);
      }
    }

118
    // For the second dim case, write also the interpolation points.
119
    if ((dim_ == 2) && (degree_ > 1)) {
Thomas Witkowski's avatar
Thomas Witkowski committed
120
121
      DOFVector< std::list< WorldVector<double> > > *interpPointCoords = (*dc_)[0]->getInterpPointCoords();
      DOFVector< std::list< WorldVector<double> > >::Iterator pointIt(interpPointCoords, USED_DOFS);
122
123
      
      for (pointIt.reset(); !pointIt.end(); ++pointIt) {
Thomas Witkowski's avatar
Thomas Witkowski committed
124
	std::list< WorldVector<double> >::iterator it2;
125
126
127
	for (it2 = pointIt->begin(); it2 != pointIt->end(); ++it2) {
	  writeCoord(file, *it2);
	}
128
129
130
131
132
      }
    }
  }


Thomas Witkowski's avatar
Thomas Witkowski committed
133
  void VtkWriter::writeVertexValues(std::ofstream &file, int componentNo)
134
  { 
135
136
    //    DOFVector<int> *interpPointInd;
    //    DOFVector<double> *values;
Thomas Witkowski's avatar
Thomas Witkowski committed
137
    //    DOFVector< std::list<WorldVector<double> > > *dofCoords;
138

139
140
    DOFVector<int> *interpPointInd = (*dc_)[componentNo]->getInterpPointInd();
    DOFVector<double> *values = (*dc_)[componentNo]->getValues();
Thomas Witkowski's avatar
Thomas Witkowski committed
141
    DOFVector< std::list<WorldVector<double> > > *dofCoords = (*dc_)[componentNo]->getDofCoords();
142
143

    /*
Thomas Witkowski's avatar
Thomas Witkowski committed
144
#ifdef _OPENMP
145
#pragma omp critical 
Thomas Witkowski's avatar
Thomas Witkowski committed
146
#endif
147
148
149
150
151
    {   
      interpPointInd = (*dc_)[componentNo]->getInterpPointInd();
      values = (*dc_)[componentNo]->getValues();
      dofCoords = (*dc_)[componentNo]->getDofCoords();
    }
152
    */
153
154
    DOFVector<int>::Iterator intPointIt(interpPointInd, USED_DOFS);
    DOFVector<double>::Iterator valueIt(values, USED_DOFS);
Thomas Witkowski's avatar
Thomas Witkowski committed
155
    DOFVector< std::list<WorldVector<double> > >::Iterator coordIt(dofCoords, USED_DOFS);
156
          
157
    // Write the values for all vertex DOFs.
158
159
160
161
162
163
    for (intPointIt.reset(), valueIt.reset(), coordIt.reset();
	 !intPointIt.end(); 
	 ++intPointIt, ++valueIt, ++coordIt) {

      if (*intPointIt == -2) {
	for (int i = 0; i < static_cast<int>(coordIt->size()); i++) {
Thomas Witkowski's avatar
Thomas Witkowski committed
164
	  file << " " << (fabs(*valueIt) < 1e-40 ? 0.0 : *valueIt) << "\n";
165
166
167
168
	}
      }
    }        

169
    // For the second dim case, write also the values of the interpolation points.
170
    if ((dim_ == 2) && (degree_ > 1)) {
Thomas Witkowski's avatar
Thomas Witkowski committed
171
      DOFVector< std::list<WorldVector<double> > >::Iterator 
172
173
174
175
176
177
178
179
	interpCoordIt((*dc_)[componentNo]->getInterpPointCoords(), USED_DOFS);
    
      for (intPointIt.reset(), valueIt.reset(), interpCoordIt.reset();
	   !intPointIt.end(); 
	   ++intPointIt, ++valueIt, ++interpCoordIt) {      
	
	if (*intPointIt >= 0) {
	  for (int i = 0; i < static_cast<int>(interpCoordIt->size()); i++) {
Thomas Witkowski's avatar
Thomas Witkowski committed
180
	    file << " " << (fabs(*valueIt) < 1e-40 ? 0.0 : *valueIt) << "\n";
181
	  }
182
183
	}
      }
184
    }    
185
186
187
  }


Thomas Witkowski's avatar
Thomas Witkowski committed
188
  void VtkWriter::writeConnectivity(std::ofstream &file) 
189
  {
190
191
    // For the second dim case, and if higher order Lagrange elements are used,
    // write the connectivity by extra functions.
192
193
194
195
196
197
198
    if ((dim_ == 2) && (degree_ == 2)) {
      writeConnectivity_dim2_degree2(file);
    } else if ((dim_ == 2) && (degree_ == 3)) {
      writeConnectivity_dim2_degree3(file);
    } else if ((dim_ == 2) && (degree_ == 4)) {
      writeConnectivity_dim2_degree4(file);   
    } else {
Thomas Witkowski's avatar
Thomas Witkowski committed
199
200
      std::list<ElementInfo> *elements = (*dc_)[0]->getElementInfos();
      std::list<ElementInfo>::iterator elementIt;
201
202
203
204
205
206
207
      int vertices = (*dc_)[0]->getMesh()->getGeo(VERTEX);
      
      for (elementIt = elements->begin(); elementIt != elements->end(); ++elementIt) {
	// for all vertices
	for (int i = 0; i < vertices; i++) {
	  file << " " << elementIt->vertexInfo[i]->outputIndex;
	}
Thomas Witkowski's avatar
Thomas Witkowski committed
208
	file << "\n";
209
210
      } 
    }
211
212
213
  }


Thomas Witkowski's avatar
Thomas Witkowski committed
214
  void VtkWriter::writeConnectivity_dim2_degree2(std::ofstream &file)
215
  {  
Thomas Witkowski's avatar
Thomas Witkowski committed
216
217
    std::list<ElementInfo> *elements = (*dc_)[0]->getElementInfos();
    std::list<ElementInfo>::iterator elementIt;
218

Thomas Witkowski's avatar
Thomas Witkowski committed
219
220
    std::vector< std::vector<int> > *interpPoints = (*dc_)[0]->getInterpPoints();
    std::vector< std::vector<int> >::iterator pointIt;
221
222
223
224
225
226
227
228
229

    int nVertices = (*dc_)[0]->getNumberVertices();

    for (pointIt = interpPoints->begin(), elementIt = elements->begin(); 
	 pointIt != interpPoints->end(); 
	 ++pointIt, ++elementIt) {
    
      file << " " << elementIt->vertexInfo[0]->outputIndex
	   << " " << (*pointIt)[1] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
230
	   << " " << (*pointIt)[2] + nVertices << "\n";
231
232
233

      file << " " << elementIt->vertexInfo[2]->outputIndex 
	   << " " << (*pointIt)[0] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
234
	   << " " << (*pointIt)[1] + nVertices << "\n";
235
236
237

      file << " " << elementIt->vertexInfo[1]->outputIndex
	   << " " << (*pointIt)[0] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
238
	   << " " << (*pointIt)[2] + nVertices << "\n";
239
240
241

      file << " " << (*pointIt)[0] + nVertices
	   << " " << (*pointIt)[1] + nVertices 
Thomas Witkowski's avatar
Thomas Witkowski committed
242
	   << " " << (*pointIt)[2] + nVertices << "\n";
243
244
245
246
247
    }

  }


Thomas Witkowski's avatar
Thomas Witkowski committed
248
  void VtkWriter::writeConnectivity_dim2_degree3(std::ofstream &file)
249
  {  
Thomas Witkowski's avatar
Thomas Witkowski committed
250
251
    std::list<ElementInfo> *elements = (*dc_)[0]->getElementInfos();
    std::list<ElementInfo>::iterator elementIt;
252

Thomas Witkowski's avatar
Thomas Witkowski committed
253
254
    std::vector< std::vector<int> > *interpPoints = (*dc_)[0]->getInterpPoints();
    std::vector< std::vector<int> >::iterator pointIt;
255
256
257
258
259
260
261
262
263

    int nVertices = (*dc_)[0]->getNumberVertices();

    for (pointIt = interpPoints->begin(), elementIt = elements->begin(); 
	 pointIt != interpPoints->end(); 
	 ++pointIt, ++elementIt) {
                          
      file << " " << elementIt->vertexInfo[0]->outputIndex
	   << " " << (*pointIt)[3] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
264
	   << " " << (*pointIt)[4] + nVertices << "\n";
265
266
267

      file << " " << (*pointIt)[4] + nVertices
	   << " " << (*pointIt)[5] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
268
	   << " " << (*pointIt)[6] + nVertices << "\n";
269
270
271

      file << " " << (*pointIt)[3] + nVertices
	   << " " << (*pointIt)[4] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
272
	   << " " << (*pointIt)[6] + nVertices << "\n";
273
274
275

      file << " " << (*pointIt)[2] + nVertices
	   << " " << (*pointIt)[3] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
276
	   << " " << (*pointIt)[6] + nVertices << "\n";
277
278
279

      file << " " << elementIt->vertexInfo[1]->outputIndex
	   << " " << (*pointIt)[0] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
280
	   << " " << (*pointIt)[5] + nVertices << "\n";
281
282
283

      file << " " << (*pointIt)[0] + nVertices
	   << " " << (*pointIt)[6] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
284
	   << " " << (*pointIt)[5] + nVertices << "\n";
285
286
287

      file << " " << (*pointIt)[0] + nVertices
	   << " " << (*pointIt)[1] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
288
	   << " " << (*pointIt)[6] + nVertices << "\n";
289
290
291

      file << " " << (*pointIt)[1] + nVertices
	   << " " << (*pointIt)[2] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
292
	   << " " << (*pointIt)[6] + nVertices << "\n";
293
294
295

      file << " " << elementIt->vertexInfo[2]->outputIndex 
	   << " " << (*pointIt)[1] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
296
	   << " " << (*pointIt)[2] + nVertices << "\n";
297
298
299
300
301
302
303
     
      
    }

  }


Thomas Witkowski's avatar
Thomas Witkowski committed
304
  void VtkWriter::writeConnectivity_dim2_degree4(std::ofstream &file)
305
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
306
307
    std::list<ElementInfo> *elements = (*dc_)[0]->getElementInfos();
    std::list<ElementInfo>::iterator elementIt;
308

Thomas Witkowski's avatar
Thomas Witkowski committed
309
310
    std::vector< std::vector<int> > *interpPoints = (*dc_)[0]->getInterpPoints();
    std::vector< std::vector<int> >::iterator pointIt;
311
312
313
314
315
316
317
318
319

    int nVertices = (*dc_)[0]->getNumberVertices();

    for (pointIt = interpPoints->begin(), elementIt = elements->begin(); 
	 pointIt != interpPoints->end(); 
	 ++pointIt, ++elementIt) {

      file << " " << elementIt->vertexInfo[0]->outputIndex 
	   << " " << (*pointIt)[5] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
320
	   << " " << (*pointIt)[6] + nVertices << "\n";
321
322
323

      file << " " << (*pointIt)[5] + nVertices
	   << " " << (*pointIt)[9] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
324
	   << " " << (*pointIt)[6] + nVertices << "\n";
325
326
327

      file << " " << (*pointIt)[6] + nVertices
	   << " " << (*pointIt)[7] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
328
	   << " " << (*pointIt)[9] + nVertices << "\n";
329
330
331

      file << " " << (*pointIt)[7] + nVertices
	   << " " << (*pointIt)[9] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
332
	   << " " << (*pointIt)[10] + nVertices << "\n";
333
334
335

      file << " " << (*pointIt)[7] + nVertices
	   << " " << (*pointIt)[8] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
336
	   << " " << (*pointIt)[10] + nVertices << "\n";
337
338
339

      file << " " << (*pointIt)[0] + nVertices
	   << " " << (*pointIt)[8] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
340
	   << " " << (*pointIt)[10] + nVertices << "\n";
341
342
343

      file << " " << elementIt->vertexInfo[1]->outputIndex 
	   << " " << (*pointIt)[0] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
344
	   << " " << (*pointIt)[8] + nVertices << "\n";
345
346
347

      file << " " << (*pointIt)[4] + nVertices
	   << " " << (*pointIt)[5] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
348
	   << " " << (*pointIt)[9] + nVertices << "\n";
349
350
351

      file << " " << (*pointIt)[4] + nVertices
	   << " " << (*pointIt)[9] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
352
	   << " " << (*pointIt)[11] + nVertices << "\n";
353
354
355

      file << " " << (*pointIt)[9] + nVertices
	   << " " << (*pointIt)[10] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
356
	   << " " << (*pointIt)[11] + nVertices << "\n";
357
358
359

      file << " " << (*pointIt)[1] + nVertices
	   << " " << (*pointIt)[10] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
360
	   << " " << (*pointIt)[11] + nVertices << "\n";
361
362
363

      file << " " << (*pointIt)[0] + nVertices
	   << " " << (*pointIt)[1] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
364
	   << " " << (*pointIt)[10] + nVertices << "\n";
365
366
367

      file << " " << (*pointIt)[3] + nVertices
	   << " " << (*pointIt)[4] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
368
	   << " " << (*pointIt)[11] + nVertices << "\n";
369
370
371

      file << " " << (*pointIt)[2] + nVertices
	   << " " << (*pointIt)[3] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
372
	   << " " << (*pointIt)[11] + nVertices << "\n";
373
374
375

      file << " " << (*pointIt)[1] + nVertices
	   << " " << (*pointIt)[2] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
376
	   << " " << (*pointIt)[11] + nVertices << "\n";
377
378
379

      file << " " << elementIt->vertexInfo[2]->outputIndex 
	   << " " << (*pointIt)[2] + nVertices
Thomas Witkowski's avatar
Thomas Witkowski committed
380
	   << " " << (*pointIt)[3] + nVertices << "\n";
381
382
383
384

    }
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
385
386
  int VtkWriter::updateAnimationFile(std::string valueFilename,
				     std::vector< std::string > *paraViewAnimationFrames,
387
388
389
390
391
				     const char *animationFilename)
  {
    size_t found = valueFilename.find_last_of("/\\");
    paraViewAnimationFrames->push_back(valueFilename.substr(found + 1));

Thomas Witkowski's avatar
Thomas Witkowski committed
392
    std::ofstream file;
393
394
    file.open(animationFilename);

Thomas Witkowski's avatar
Thomas Witkowski committed
395
396
397
    file << "<?xml version=\"1.0\"?>\n";
    file << "<VTKFile type=\"Collection\" version=\"0.1\" >"  << "\n";
    file << "<Collection>\n";
398
399

    int counter = 0;
Thomas Witkowski's avatar
Thomas Witkowski committed
400
    std::vector< std::string >::iterator it;
401
402
403
404
    for (it = paraViewAnimationFrames->begin(); 
	 it < paraViewAnimationFrames->end(); 
	 ++it, counter++) {
      file << "<DataSet timestep=\"" << counter
Thomas Witkowski's avatar
Thomas Witkowski committed
405
	   << "\" part=\"0\" file=\"" << (*it) << "\"/>\n";      
406
407
    }

Thomas Witkowski's avatar
Thomas Witkowski committed
408
409
    file << "</Collection>\n";
    file << "</VTKFile>\n";
410
411
412
413
414

    file.close();
    
    return 0;
  }
415
416


417
}