Dune Core Modules (2.4.1)

vtkwriter.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3
4#ifndef DUNE_VTKWRITER_HH
5#define DUNE_VTKWRITER_HH
6
7#include <cstring>
8#include <iostream>
9#include <string>
10#include <fstream>
11#include <sstream>
12#include <iomanip>
13
14#include <vector>
15#include <list>
16
20#include <dune/common/std/memory.hh>
21#include <dune/common/indent.hh>
23#include <dune/common/path.hh>
25#include <dune/geometry/referenceelements.hh>
27#include <dune/grid/common/gridenums.hh>
31#include <dune/grid/io/file/vtk/pvtuwriter.hh>
32#include <dune/grid/io/file/vtk/streams.hh>
33#include <dune/grid/io/file/vtk/vtuwriter.hh>
34
48namespace Dune
49{
50
51 namespace detail {
52
53 template<typename F, typename = int>
54 struct _has_local_context
55 : public std::false_type
56 {};
57
58 template<typename T>
59 struct _has_local_context<T,typename std::enable_if<(sizeof(static_cast<T*>(nullptr)->localContext()) > 0),int>::type>
60 : public std::true_type
61 {};
62
63 }
64
65 namespace VTKWriteTypeTraits {
66 template<typename T>
67 struct IsLocalFunction
68 {
69 };
70 }
71
72 // Forward-declaration here, so the class can be friend of VTKWriter
73 template <class GridView>
74 class VTKSequenceWriterBase;
75 template <class GridView>
76 class VTKSequenceWriter;
77
86 template< class GridView >
87 class VTKWriter {
88
89 // VTKSequenceWriterBase needs getSerialPieceName
90 // and getParallelHeaderName
91 friend class VTKSequenceWriterBase<GridView>;
92 // VTKSequenceWriter needs the grid view, to get the MPI size and rank
93 friend class VTKSequenceWriter<GridView>;
94
95 // extract types
96 typedef typename GridView::Grid Grid;
97 typedef typename GridView::ctype DT;
98 enum { n = GridView::dimension };
99 enum { w = GridView::dimensionworld };
100
101 typedef typename GridView::template Codim< 0 >::Entity Cell;
102 typedef typename GridView::template Codim< n >::Entity Vertex;
103 typedef Cell Entity;
104
105 typedef typename GridView::IndexSet IndexSet;
106
107 static const PartitionIteratorType VTK_Partition = InteriorBorder_Partition;
108 //static const PartitionIteratorType VTK_Partition = All_Partition;
109
110 typedef typename GridView::template Codim< 0 >
111 ::template Partition< VTK_Partition >::Iterator
112 GridCellIterator;
113 typedef typename GridView::template Codim< n >
114 ::template Partition< VTK_Partition >::Iterator
115 GridVertexIterator;
116
117 typedef typename GridCellIterator::Reference EntityReference;
118
119 typedef typename GridView::template Codim< 0 >
120 ::Entity::Geometry::LocalCoordinate Coordinate;
121
123
124 // return true if entity should be skipped in Vertex and Corner iterator
125 static bool skipEntity( const PartitionType entityType )
126 {
127 switch( VTK_Partition )
128 {
129 // for All_Partition no entity has to be skipped
130 case All_Partition: return false;
131 case InteriorBorder_Partition: return ( entityType != InteriorEntity );
132 default: DUNE_THROW(NotImplemented,"Add check for this partition type");
133 }
134 return false ;
135 }
136
137 public:
138
140 typedef shared_ptr< const VTKFunction > VTKFunctionPtr;
141
142 protected:
143
145
149 {
150
151 public:
152
154
157 {
158
160 virtual void bind(const Entity& e) = 0;
161
163 virtual void unbind() = 0;
164
166
169 virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const = 0;
170
171 virtual ~FunctionWrapperBase()
172 {}
173
174 };
175
177 template<typename F>
179 : public FunctionWrapperBase
180 {
181 typedef typename std::decay<F>::type Function;
182
183 template<typename F_>
184 FunctionWrapper(F_&& f)
185 : _f(std::forward<F_>(f))
186 {}
187
188 virtual void bind(const Entity& e)
189 {
190 _f.bind(e);
191 }
192
193 virtual void unbind()
194 {
195 _f.unbind();
196 }
197
198 virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
199 {
200 auto r = _f(pos);
201 // we need to do different things here depending on whether r supports indexing into it or not.
202 do_write(w,r,count,is_indexable<decltype(r)>());
203 }
204
205 private:
206
207 template<typename R>
208 void do_write(Writer& w, const R& r, std::size_t count, std::true_type) const
209 {
210 for (std::size_t i = 0; i < count; ++i)
211 w.write(r[i]);
212 }
213
214 template<typename R>
215 void do_write(Writer& w, const R& r, std::size_t count, std::false_type) const
216 {
217 assert(count == 1);
218 w.write(r);
219 }
220
221 Function _f;
222 };
223
226 : public FunctionWrapperBase
227 {
228 VTKFunctionWrapper(const VTKFunctionPtr& f)
229 : _f(f)
230 , _entity(nullptr)
231 {}
232
233 virtual void bind(const Entity& e)
234 {
235 _entity = &e;
236 }
237
238 virtual void unbind()
239 {
240 _entity = nullptr;
241 }
242
243 virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
244 {
245 for (std::size_t i = 0; i < count; ++i)
246 w.write(_f->evaluate(i,*_entity,pos));
247 }
248
249 private:
250
251 VTKFunctionPtr _f;
252 const Entity* _entity;
253
254 };
255
257 template<typename F>
259 typename std::enable_if<detail::_has_local_context<F>::value,int>::type dummy = 0)
260 : _f(Dune::Std::make_unique<FunctionWrapper<F> >(std::forward<F>(f)))
261 , _fieldInfo(fieldInfo)
262 {}
263
265 template<typename F>
267 typename std::enable_if<not detail::_has_local_context<F>::value,int>::type dummy = 0)
268 : _f(Dune::Std::make_unique< FunctionWrapper<
269 typename std::decay<decltype(localFunction(std::forward<F>(f)))>::type
270 > >(localFunction(std::forward<F>(f))))
271 , _fieldInfo(fieldInfo)
272 {}
273
275 explicit VTKLocalFunction (const VTKFunctionPtr& vtkFunctionPtr)
276 : _f(Dune::Std::make_unique<VTKFunctionWrapper>(vtkFunctionPtr))
277 , _fieldInfo(
278 vtkFunctionPtr->name(),
279 vtkFunctionPtr->ncomps() > 1 ? VTK::FieldInfo::Type::vector : VTK::FieldInfo::Type::scalar,
280 vtkFunctionPtr->ncomps()
281 )
282 {}
283
285 std::string name() const
286 {
287 return fieldInfo().name();
288 }
289
292 {
293 return _fieldInfo;
294 }
295
297 void bind(const Entity& e) const
298 {
299 _f->bind(e);
300 }
301
303 void unbind() const
304 {
305 _f->unbind();
306 }
307
309 void write(const Coordinate& pos, Writer& w) const
310 {
311 _f->write(pos,w,fieldInfo().size());
312 }
313
314 std::shared_ptr<FunctionWrapperBase> _f;
315 VTK::FieldInfo _fieldInfo;
316
317 };
318
319 typedef typename std::list<VTKLocalFunction>::const_iterator FunctionIterator;
320
322
327 class CellIterator : public GridCellIterator
328 {
329 public:
331 CellIterator(const GridCellIterator & x) : GridCellIterator(x) {}
335 {
336 return ReferenceElements<DT,n>::general((*this)->type()).position(0,0);
337 }
338 };
339
340 CellIterator cellBegin() const
341 {
342 return gridView_.template begin< 0, VTK_Partition >();
343 }
344
345 CellIterator cellEnd() const
346 {
347 return gridView_.template end< 0, VTK_Partition >();
348 }
349
351
366 public ForwardIteratorFacade<VertexIterator, const Entity, EntityReference, int>
367 {
368 GridCellIterator git;
369 GridCellIterator gend;
370 VTK::DataMode datamode;
371 // Index of the currently visited corner within the current element.
372 // NOTE: this is in Dune-numbering, in contrast to CornerIterator.
373 int cornerIndexDune;
374 const VertexMapper & vertexmapper;
375 std::vector<bool> visited;
376 // in conforming mode, for each vertex id (as obtained by vertexmapper)
377 // hold its number in the iteration order (VertexIterator)
378 int offset;
379
380 // hide operator ->
381 void operator->();
382 protected:
383 void basicIncrement ()
384 {
385 if( git == gend )
386 return;
387 ++cornerIndexDune;
388 const int numCorners = git->subEntities(n);
389 if( cornerIndexDune == numCorners )
390 {
391 offset += numCorners;
392 cornerIndexDune = 0;
393
394 ++git;
395 while( (git != gend) && skipEntity( git->partitionType() ) )
396 ++git;
397 }
398 }
399 public:
400 VertexIterator(const GridCellIterator & x,
401 const GridCellIterator & end,
402 const VTK::DataMode & dm,
403 const VertexMapper & vm) :
404 git(x), gend(end), datamode(dm), cornerIndexDune(0),
405 vertexmapper(vm), visited(vm.size(), false),
406 offset(0)
407 {
408 if (datamode == VTK::conforming && git != gend)
409 visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
410 }
411 void increment ()
412 {
413 switch (datamode)
414 {
415 case VTK::conforming :
416 while(visited[vertexmapper.subIndex(*git,cornerIndexDune,n)])
417 {
418 basicIncrement();
419 if (git == gend) return;
420 }
421 visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
422 break;
423 case VTK::nonconforming :
424 basicIncrement();
425 break;
426 }
427 }
428 bool equals (const VertexIterator & cit) const
429 {
430 return git == cit.git
431 && cornerIndexDune == cit.cornerIndexDune
432 && datamode == cit.datamode;
433 }
434 EntityReference dereference() const
435 {
436 return *git;
437 }
439 int localindex () const
440 {
441 return cornerIndexDune;
442 }
445 {
446 return ReferenceElements<DT,n>::general(git->type())
447 .position(cornerIndexDune,n);
448 }
449 };
450
451 VertexIterator vertexBegin () const
452 {
453 return VertexIterator( gridView_.template begin< 0, VTK_Partition >(),
454 gridView_.template end< 0, VTK_Partition >(),
455 datamode, *vertexmapper );
456 }
457
458 VertexIterator vertexEnd () const
459 {
460 return VertexIterator( gridView_.template end< 0, VTK_Partition >(),
461 gridView_.template end< 0, VTK_Partition >(),
462 datamode, *vertexmapper );
463 }
464
466
481 public ForwardIteratorFacade<CornerIterator, const Entity, EntityReference, int>
482 {
483 GridCellIterator git;
484 GridCellIterator gend;
485 VTK::DataMode datamode;
486 // Index of the currently visited corner within the current element.
487 // NOTE: this is in VTK-numbering, in contrast to VertexIterator.
488 int cornerIndexVTK;
489 const VertexMapper & vertexmapper;
490 // in conforming mode, for each vertex id (as obtained by vertexmapper)
491 // hold its number in the iteration order of VertexIterator (*not*
492 // CornerIterator)
493 const std::vector<int> & number;
494 // holds the number of corners of all the elements we have seen so far,
495 // excluding the current element
496 int offset;
497
498 // hide operator ->
499 void operator->();
500 public:
501 CornerIterator(const GridCellIterator & x,
502 const GridCellIterator & end,
503 const VTK::DataMode & dm,
504 const VertexMapper & vm,
505 const std::vector<int> & num) :
506 git(x), gend(end), datamode(dm), cornerIndexVTK(0),
507 vertexmapper(vm),
508 number(num), offset(0) {}
509 void increment ()
510 {
511 if( git == gend )
512 return;
513 ++cornerIndexVTK;
514 const int numCorners = git->subEntities(n);
515 if( cornerIndexVTK == numCorners )
516 {
517 offset += numCorners;
518 cornerIndexVTK = 0;
519
520 ++git;
521 while( (git != gend) && skipEntity( git->partitionType() ) )
522 ++git;
523 }
524 }
525 bool equals (const CornerIterator & cit) const
526 {
527 return git == cit.git
528 && cornerIndexVTK == cit.cornerIndexVTK
529 && datamode == cit.datamode;
530 }
531 EntityReference dereference() const
532 {
533 return *git;
534 }
536
540 int id () const
541 {
542 switch (datamode)
543 {
544 case VTK::conforming :
545 return
546 number[vertexmapper.subIndex(*git,VTK::renumber(*git,cornerIndexVTK),
547 n)];
548 case VTK::nonconforming :
549 return offset + VTK::renumber(*git,cornerIndexVTK);
550 default :
551 DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
552 }
553 }
554 };
555
556 CornerIterator cornerBegin () const
557 {
558 return CornerIterator( gridView_.template begin< 0, VTK_Partition >(),
559 gridView_.template end< 0, VTK_Partition >(),
560 datamode, *vertexmapper, number );
561 }
562
563 CornerIterator cornerEnd () const
564 {
565 return CornerIterator( gridView_.template end< 0, VTK_Partition >(),
566 gridView_.template end< 0, VTK_Partition >(),
567 datamode, *vertexmapper, number );
568 }
569
570 public:
578 explicit VTKWriter ( const GridView &gridView,
579 VTK::DataMode dm = VTK::conforming )
580 : gridView_( gridView ),
581 datamode( dm )
582 { }
583
588 void addCellData (const VTKFunctionPtr & p)
589 {
590 celldata.push_back(VTKLocalFunction(p));
591 }
592
593 template<typename F>
594 void addCellData(F&& f, VTK::FieldInfo vtkFieldInfo)
595 {
596 celldata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
597 }
598
604 void addCellData (VTKFunction* p) DUNE_DEPRECATED_MSG("Don't pass raw pointers, use the version with shared_ptr")
605 {
606 celldata.push_back(VTKLocalFunction(VTKFunctionPtr(p)));
607 }
608
624 template<class V>
625 void addCellData (const V& v, const std::string &name, int ncomps = 1)
626 {
628 for (int c=0; c<ncomps; ++c) {
629 std::stringstream compName;
630 compName << name;
631 if (ncomps>1)
632 compName << "[" << c << "]";
633 VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
634 addCellData(VTKFunctionPtr(p));
635 }
636 }
637
643 void addVertexData (VTKFunction* p) DUNE_DEPRECATED_MSG("Don't pass raw pointers, use the version with shared_ptr")
644 {
645 vertexdata.push_back(VTKLocalFunction(VTKFunctionPtr(p)));
646 }
647
652 void addVertexData (const VTKFunctionPtr & p)
653 {
654 vertexdata.push_back(VTKLocalFunction(p));
655 }
656
657 template<typename F>
658 void addVertexData(F&& f, VTK::FieldInfo vtkFieldInfo)
659 {
660 vertexdata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
661 }
662
663
679 template<class V>
680 void addVertexData (const V& v, const std::string &name, int ncomps=1)
681 {
683 for (int c=0; c<ncomps; ++c) {
684 std::stringstream compName;
685 compName << name;
686 if (ncomps>1)
687 compName << "[" << c << "]";
688 VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
689 addVertexData(VTKFunctionPtr(p));
690 }
691 }
692
694 void clear ()
695 {
696 celldata.clear();
697 vertexdata.clear();
698 }
699
701 virtual ~VTKWriter ()
702 {
703 this->clear();
704 }
705
717 std::string write ( const std::string &name,
718 VTK::OutputType type = VTK::ascii )
719 {
720 return write( name, type, gridView_.comm().rank(), gridView_.comm().size() );
721 }
722
749 std::string pwrite ( const std::string & name, const std::string & path, const std::string & extendpath,
750 VTK::OutputType type = VTK::ascii )
751 {
752 return pwrite( name, path, extendpath, type, gridView_.comm().rank(), gridView_.comm().size() );
753 }
754
755 protected:
757
768 std::string getParallelPieceName(const std::string& name,
769 const std::string& path,
770 int commRank, int commSize) const
771 {
772 std::ostringstream s;
773 if(path.size() > 0) {
774 s << path;
775 if(path[path.size()-1] != '/')
776 s << '/';
777 }
778 s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
779 s << 'p' << std::setw(4) << std::setfill('0') << commRank << '-';
780 s << name;
781 if(GridView::dimension > 1)
782 s << ".vtu";
783 else
784 s << ".vtp";
785 return s.str();
786 }
787
789
799 std::string getParallelHeaderName(const std::string& name,
800 const std::string& path,
801 int commSize) const
802 {
803 std::ostringstream s;
804 if(path.size() > 0) {
805 s << path;
806 if(path[path.size()-1] != '/')
807 s << '/';
808 }
809 s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
810 s << name;
811 if(GridView::dimension > 1)
812 s << ".pvtu";
813 else
814 s << ".pvtp";
815 return s.str();
816 }
817
819
831 std::string getSerialPieceName(const std::string& name,
832 const std::string& path) const
833 {
834 static const std::string extension =
835 GridView::dimension == 1 ? ".vtp" : ".vtu";
836
837 return concatPaths(path, name+extension);
838 }
839
855 std::string write ( const std::string &name,
856 VTK::OutputType type,
857 const int commRank,
858 const int commSize )
859 {
860 // in the parallel case, just use pwrite, it has all the necessary
861 // stuff, so we don't need to reimplement it here.
862 if(commSize > 1)
863 return pwrite(name, "", "", type, commRank, commSize);
864
865 // make data mode visible to private functions
866 outputtype = type;
867
868 // generate filename for process data
869 std::string pieceName = getSerialPieceName(name, "");
870
871 // write process data
872 std::ofstream file;
873 file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
874 std::ios_base::eofbit);
875 // check if file can be opened
876 try {
877 file.open( pieceName.c_str(), std::ios::binary );
878 }
879 catch(...) {
880 std::cerr << "Filename: " << pieceName << " could not be opened" << std::endl;
881 throw;
882 }
883 if (! file.is_open())
884 DUNE_THROW(IOError, "Could not write to piece file " << pieceName);
885 writeDataFile( file );
886 file.close();
887
888 return pieceName;
889 }
890
892
915 std::string pwrite(const std::string& name, const std::string& path,
916 const std::string& extendpath,
917 VTK::OutputType ot, const int commRank,
918 const int commSize )
919 {
920 // make data mode visible to private functions
921 outputtype=ot;
922
923 // do some magic because paraview can only cope with relative pathes to piece files
924 std::ofstream file;
925 file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
926 std::ios_base::eofbit);
927 std::string piecepath = concatPaths(path, extendpath);
928 std::string relpiecepath = relativePath(path, piecepath);
929
930 // write this processes .vtu/.vtp piece file
931 std::string fullname = getParallelPieceName(name, piecepath, commRank,
932 commSize);
933 // check if file can be opened
934 try {
935 file.open(fullname.c_str(),std::ios::binary);
936 }
937 catch(...) {
938 std::cerr << "Filename: " << fullname << " could not be opened" << std::endl;
939 throw;
940 }
941 if (! file.is_open())
942 DUNE_THROW(IOError, "Could not write to piecefile file " << fullname);
943 writeDataFile(file);
944 file.close();
945 gridView_.comm().barrier();
946
947 // if we are rank 0, write .pvtu/.pvtp parallel header
948 fullname = getParallelHeaderName(name, path, commSize);
949 if( commRank ==0 )
950 {
951 file.open(fullname.c_str());
952 if (! file.is_open())
953 DUNE_THROW(IOError, "Could not write to parallel file " << fullname);
954 writeParallelHeader(file,name,relpiecepath, commSize );
955 file.close();
956 }
957 gridView_.comm().barrier();
958 return fullname;
959 }
960
961 private:
963
980 void writeParallelHeader(std::ostream& s, const std::string& piecename,
981 const std::string& piecepath, const int commSize)
982 {
983 VTK::FileType fileType =
984 (n == 1) ? VTK::polyData : VTK::unstructuredGrid;
985
986 VTK::PVTUWriter writer(s, fileType);
987
988 writer.beginMain();
989
990 // PPointData
991 {
992 std::string scalars, vectors;
993 std::tie(scalars,vectors) = getDataNames(vertexdata);
994 writer.beginPointData(scalars, vectors);
995 }
996 for (auto it = vertexdata.begin(),
997 end = vertexdata.end();
998 it != end;
999 ++it)
1000 {
1001 unsigned writecomps = it->fieldInfo().size();
1002 if(writecomps == 2) writecomps = 3;
1003 writer.addArray<float>(it->name(), writecomps);
1004 }
1005 writer.endPointData();
1006
1007 // PCellData
1008 {
1009 std::string scalars, vectors;
1010 std::tie(scalars,vectors) = getDataNames(celldata);
1011 writer.beginCellData(scalars, vectors);
1012 }
1013 for (auto it = celldata.begin(),
1014 end = celldata.end();
1015 it != end;
1016 ++it)
1017 {
1018 unsigned writecomps = it->fieldInfo().size();
1019 if(writecomps == 2) writecomps = 3;
1020 writer.addArray<float>(it->name(), writecomps);
1021 }
1022 writer.endCellData();
1023
1024 // PPoints
1025 writer.beginPoints();
1026 writer.addArray<float>("Coordinates", 3);
1027 writer.endPoints();
1028
1029 // Pieces
1030 for( int i = 0; i < commSize; ++i )
1031 {
1032 const std::string& fullname = getParallelPieceName(piecename,
1033 piecepath, i,
1034 commSize);
1035 writer.addPiece(fullname);
1036 }
1037
1038 writer.endMain();
1039 }
1040
1042 void writeDataFile (std::ostream& s)
1043 {
1044 VTK::FileType fileType =
1045 (n == 1) ? VTK::polyData : VTK::unstructuredGrid;
1046
1047 VTK::VTUWriter writer(s, outputtype, fileType);
1048
1049 // Grid characteristics
1050 vertexmapper = new VertexMapper( gridView_ );
1051 if (datamode == VTK::conforming)
1052 {
1053 number.resize(vertexmapper->size());
1054 for (std::vector<int>::size_type i=0; i<number.size(); i++) number[i] = -1;
1055 }
1056 countEntities(nvertices, ncells, ncorners);
1057
1058 writer.beginMain(ncells, nvertices);
1059 writeAllData(writer);
1060 writer.endMain();
1061
1062 // write appended binary data section
1063 if(writer.beginAppended())
1064 writeAllData(writer);
1065 writer.endAppended();
1066
1067 delete vertexmapper; number.clear();
1068 }
1069
1070 void writeAllData(VTK::VTUWriter& writer) {
1071 // PointData
1072 writeVertexData(writer);
1073
1074 // CellData
1075 writeCellData(writer);
1076
1077 // Points
1078 writeGridPoints(writer);
1079
1080 // Cells
1081 writeGridCells(writer);
1082 }
1083
1084 protected:
1085 std::string getFormatString() const
1086 {
1087 if (outputtype==VTK::ascii)
1088 return "ascii";
1089 if (outputtype==VTK::base64)
1090 return "binary";
1091 if (outputtype==VTK::appendedraw)
1092 return "appended";
1093 if (outputtype==VTK::appendedbase64)
1094 return "appended";
1095 DUNE_THROW(IOError, "VTKWriter: unsupported OutputType" << outputtype);
1096 }
1097
1098 std::string getTypeString() const
1099 {
1100 if (n==1)
1101 return "PolyData";
1102 else
1103 return "UnstructuredGrid";
1104 }
1105
1107 virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
1108 {
1109 nvertices = 0;
1110 ncells = 0;
1111 ncorners = 0;
1112 for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1113 {
1114 ncells++;
1115 // because of the use of vertexmapper->map(), this iteration must be
1116 // in the order of Dune's numbering.
1117 const int subEntities = it->subEntities(n);
1118 for (int i=0; i<subEntities; ++i)
1119 {
1120 ncorners++;
1121 if (datamode == VTK::conforming)
1122 {
1123 int alpha = vertexmapper->subIndex(*it,i,n);
1124 if (number[alpha]<0)
1125 number[alpha] = nvertices++;
1126 }
1127 else
1128 {
1129 nvertices++;
1130 }
1131 }
1132 }
1133 }
1134
1135 template<typename T>
1136 std::tuple<std::string,std::string> getDataNames(const T& data) const
1137 {
1138 std::string scalars = "";
1139 for (auto it = data.begin(),
1140 end = data.end();
1141 it != end;
1142 ++it)
1143 if (it->fieldInfo().type() == VTK::FieldInfo::Type::scalar)
1144 {
1145 scalars = it->name();
1146 break;
1147 }
1148
1149 std::string vectors = "";
1150 for (auto it = data.begin(),
1151 end = data.end();
1152 it != end;
1153 ++it)
1154 if (it->fieldInfo().type() == VTK::FieldInfo::Type::vector)
1155 {
1156 vectors = it->name();
1157 break;
1158 }
1159 return std::make_tuple(scalars,vectors);
1160 }
1161
1162 template<typename Data, typename Iterator>
1163 void writeData(VTK::VTUWriter& writer, const Data& data, const Iterator begin, const Iterator end, int nentries)
1164 {
1165 for (auto it = data.begin(),
1166 iend = data.end();
1167 it != iend;
1168 ++it)
1169 {
1170 const auto& f = *it;
1171 VTK::FieldInfo fieldInfo = f.fieldInfo();
1172 std::size_t writecomps = fieldInfo.size();
1173 switch (fieldInfo.type())
1174 {
1176 break;
1178 // vtk file format: a vector data always should have 3 comps (with
1179 // 3rd comp = 0 in 2D case)
1180 if (writecomps > 3)
1181 DUNE_THROW(IOError,"Cannot write VTK vectors with more than 3 components (components was " << writecomps << ")");
1182 writecomps = 3;
1183 break;
1185 DUNE_THROW(NotImplemented,"VTK output for tensors not implemented yet");
1186 }
1187 shared_ptr<VTK::DataArrayWriter<float> > p
1188 (writer.makeArrayWriter<float>(f.name(), writecomps, nentries));
1189 if(!p->writeIsNoop())
1190 for (Iterator eit = begin; eit!=end; ++eit)
1191 {
1192 const Entity & e = *eit;
1193 f.bind(e);
1194 f.write(eit.position(),*p);
1195 f.unbind();
1196 // vtk file format: a vector data always should have 3 comps
1197 // (with 3rd comp = 0 in 2D case)
1198 for (std::size_t j=fieldInfo.size(); j < writecomps; ++j)
1199 p->write(0.0);
1200 }
1201 }
1202 }
1203
1205 virtual void writeCellData(VTK::VTUWriter& writer)
1206 {
1207 if(celldata.size() == 0)
1208 return;
1209
1210 std::string scalars, vectors;
1211 std::tie(scalars,vectors) = getDataNames(celldata);
1212
1213 writer.beginCellData(scalars, vectors);
1214 writeData(writer,celldata,cellBegin(),cellEnd(),ncells);
1215 writer.endCellData();
1216 }
1217
1219 virtual void writeVertexData(VTK::VTUWriter& writer)
1220 {
1221 if(vertexdata.size() == 0)
1222 return;
1223
1224 std::string scalars, vectors;
1225 std::tie(scalars,vectors) = getDataNames(vertexdata);
1226
1227 writer.beginPointData(scalars, vectors);
1228 writeData(writer,vertexdata,vertexBegin(),vertexEnd(),nvertices);
1229 writer.endPointData();
1230 }
1231
1233 virtual void writeGridPoints(VTK::VTUWriter& writer)
1234 {
1235 writer.beginPoints();
1236
1237 shared_ptr<VTK::DataArrayWriter<float> > p
1238 (writer.makeArrayWriter<float>("Coordinates", 3, nvertices));
1239 if(!p->writeIsNoop()) {
1240 VertexIterator vEnd = vertexEnd();
1241 for (VertexIterator vit=vertexBegin(); vit!=vEnd; ++vit)
1242 {
1243 int dimw=w;
1244 for (int j=0; j<std::min(dimw,3); j++)
1245 p->write((*vit).geometry().corner(vit.localindex())[j]);
1246 for (int j=std::min(dimw,3); j<3; j++)
1247 p->write(0.0);
1248 }
1249 }
1250 // free the VTK::DataArrayWriter before touching the stream
1251 p.reset();
1252
1253 writer.endPoints();
1254 }
1255
1257 virtual void writeGridCells(VTK::VTUWriter& writer)
1258 {
1259 writer.beginCells();
1260
1261 // connectivity
1262 {
1263 shared_ptr<VTK::DataArrayWriter<int> > p1
1264 (writer.makeArrayWriter<int>("connectivity", 1, ncorners));
1265 if(!p1->writeIsNoop())
1266 for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
1267 p1->write(it.id());
1268 }
1269
1270 // offsets
1271 {
1272 shared_ptr<VTK::DataArrayWriter<int> > p2
1273 (writer.makeArrayWriter<int>("offsets", 1, ncells));
1274 if(!p2->writeIsNoop()) {
1275 int offset = 0;
1276 for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1277 {
1278 offset += it->subEntities(n);
1279 p2->write(offset);
1280 }
1281 }
1282 }
1283
1284 // types
1285 if (n>1)
1286 {
1287 shared_ptr<VTK::DataArrayWriter<unsigned char> > p3
1288 (writer.makeArrayWriter<unsigned char>("types", 1, ncells));
1289 if(!p3->writeIsNoop())
1290 for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1291 {
1292 int vtktype = VTK::geometryType(it->type());
1293 p3->write(vtktype);
1294 }
1295 }
1296
1297 writer.endCells();
1298 }
1299
1300 protected:
1301 // the list of registered functions
1302 std::list<VTKLocalFunction> celldata;
1303 std::list<VTKLocalFunction> vertexdata;
1304
1305 // the grid
1306 GridView gridView_;
1307
1308 // temporary grid information
1309 int ncells;
1310 int nvertices;
1311 int ncorners;
1312 private:
1313 VertexMapper* vertexmapper;
1314 // in conforming mode, for each vertex id (as obtained by vertexmapper)
1315 // hold its number in the iteration order (VertexIterator)
1316 std::vector<int> number;
1317 VTK::DataMode datamode;
1318 protected:
1319 VTK::OutputType outputtype;
1320 };
1321
1322}
1323
1324#endif
vector space out of a tensor product of fields.
Definition: fvector.hh:94
Base class for stl conformant forward iterators.
Definition: iteratorfacades.hh:141
Base class template for function classes.
Definition: function.hh:28
Grid view abstract base class.
Definition: gridview.hh:59
Default exception class for I/O errors.
Definition: exceptions.hh:256
Implementation class for a multiple codim and multiple geometry type mapper.
Definition: mcmgmapper.hh:104
Default exception for dummy implementations.
Definition: exceptions.hh:288
Take a vector and interpret it as cell data for the VTKWriter.
Definition: function.hh:90
Take a vector and interpret it as point data for the VTKWriter.
Definition: function.hh:188
A base class for grid functions with any return type and dimension.
Definition: function.hh:39
Base class to write pvd-files which contains a list of all collected vtk-files.
Definition: vtksequencewriterbase.hh:33
Writer for the ouput of grid functions in the vtk format.
Definition: vtksequencewriter.hh:25
Iterator over the grids elements.
Definition: vtkwriter.hh:328
CellIterator(const GridCellIterator &x)
construct a CellIterator from the gridview's Iterator.
Definition: vtkwriter.hh:331
const FieldVector< DT, n > position() const
Definition: vtkwriter.hh:334
Iterate over the elements' corners.
Definition: vtkwriter.hh:482
int id() const
Process-local consecutive zero-starting vertex id.
Definition: vtkwriter.hh:540
Type erasure wrapper for VTK data sets.
Definition: vtkwriter.hh:149
VTKLocalFunction(const VTKFunctionPtr &vtkFunctionPtr)
Construct a VTKLocalFunction for a legacy VTKFunction.
Definition: vtkwriter.hh:275
void unbind() const
Unbind the data set from the currently bound entity.
Definition: vtkwriter.hh:303
std::string name() const
Returns the name of the data set.
Definition: vtkwriter.hh:285
const VTK::FieldInfo & fieldInfo() const
Returns the VTK::FieldInfo for the data set.
Definition: vtkwriter.hh:291
void bind(const Entity &e) const
Bind the data set to grid entity e.
Definition: vtkwriter.hh:297
VTKLocalFunction(F &&f, VTK::FieldInfo fieldInfo, typename std::enable_if< detail::_has_local_context< F >::value, int >::type dummy=0)
Construct a VTKLocalFunction for a dune-functions style LocalFunction.
Definition: vtkwriter.hh:258
VTKLocalFunction(F &&f, VTK::FieldInfo fieldInfo, typename std::enable_if< not detail::_has_local_context< F >::value, int >::type dummy=0)
Construct a VTKLocalFunction for a dune-functions style Function.
Definition: vtkwriter.hh:266
void write(const Coordinate &pos, Writer &w) const
Write the value of the data set at local coordinate pos to the writer w.
Definition: vtkwriter.hh:309
Iterate over the grid's vertices.
Definition: vtkwriter.hh:367
const FieldVector< DT, n > & position() const
position of vertex inside the entity
Definition: vtkwriter.hh:444
int localindex() const
index of vertex within the entity, in Dune-numbering
Definition: vtkwriter.hh:439
Writer for the ouput of grid functions in the vtk format.
Definition: vtkwriter.hh:87
void addVertexData(const VTKFunctionPtr &p)
Add a grid function that lives on the vertices of the grid to the visualization.
Definition: vtkwriter.hh:652
void clear()
clear list of registered functions
Definition: vtkwriter.hh:694
std::string write(const std::string &name, VTK::OutputType type=VTK::ascii)
write output (interface might change later)
Definition: vtkwriter.hh:717
std::string getParallelHeaderName(const std::string &name, const std::string &path, int commSize) const
return name of a parallel header file
Definition: vtkwriter.hh:799
std::string getSerialPieceName(const std::string &name, const std::string &path) const
return name of a serial piece file
Definition: vtkwriter.hh:831
void addCellData(const VTKFunctionPtr &p)
Add a grid function that lives on the cells of the grid to the visualization.
Definition: vtkwriter.hh:588
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: vtkwriter.hh:1205
VTKWriter(const GridView &gridView, VTK::DataMode dm=VTK::conforming)
Construct a VTKWriter working on a specific GridView.
Definition: vtkwriter.hh:578
std::string getParallelPieceName(const std::string &name, const std::string &path, int commRank, int commSize) const
return name of a parallel piece file
Definition: vtkwriter.hh:768
virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
count the vertices, cells and corners
Definition: vtkwriter.hh:1107
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: vtkwriter.hh:1257
void addVertexData(const V &v, const std::string &name, int ncomps=1)
Add a grid function (represented by container) that lives on the vertices of the grid to the visualiz...
Definition: vtkwriter.hh:680
void addCellData(VTKFunction *p) DUNE_DEPRECATED_MSG("Don't pass raw pointers
Add a grid function that lives on the cells of the grid to the visualization.
std::string write(const std::string &name, VTK::OutputType type, const int commRank, const int commSize)
write output (interface might change later)
Definition: vtkwriter.hh:855
void addCellData(const V &v, const std::string &name, int ncomps=1)
Add a grid function (represented by container) that lives on the cells of the grid to the visualizati...
Definition: vtkwriter.hh:625
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: vtkwriter.hh:1233
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: vtkwriter.hh:1219
virtual ~VTKWriter()
destructor
Definition: vtkwriter.hh:701
void addVertexData(VTKFunction *p) DUNE_DEPRECATED_MSG("Don't pass raw pointers
Add a grid function that lives on the vertices of the grid to the visualization.
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath, VTK::OutputType ot, const int commRank, const int commSize)
write output; interface might change later
Definition: vtkwriter.hh:915
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath, VTK::OutputType type=VTK::ascii)
write output (interface might change later)
Definition: vtkwriter.hh:749
virtual void write(T data)=0
write one data element
Descriptor struct for VTK fields.
Definition: common.hh:307
@ tensor
tensor field (always 3x3)
@ vector
vector-valued field (always 3D, will be padded if necessary)
std::string name() const
The name of the data field.
Definition: common.hh:330
Dump a .vtu/.vtp files contents to a stream.
Definition: pvtuwriter.hh:60
Dump a .vtu/.vtp files contents to a stream.
Definition: vtuwriter.hh:96
void endCellData()
finish CellData section
Definition: vtuwriter.hh:218
void beginCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:272
DataArrayWriter< T > * makeArrayWriter(const std::string &name, unsigned ncomps, unsigned nitems)
aquire a DataArrayWriter
Definition: vtuwriter.hh:379
void endPointData()
finish PointData section
Definition: vtuwriter.hh:180
void beginCellData(const std::string &scalars="", const std::string &vectors="")
start CellData section
Definition: vtuwriter.hh:203
void beginPointData(const std::string &scalars="", const std::string &vectors="")
start PointData section
Definition: vtuwriter.hh:165
void endPoints()
finish section for the point coordinates
Definition: vtuwriter.hh:247
void endCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:283
void beginPoints()
start section for the point coordinates
Definition: vtuwriter.hh:236
Common stuff for the VTKWriter.
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:40
FileType
which type of VTK file to write
Definition: common.hh:290
DataMode
Whether to produce conforming or non-conforming output.
Definition: common.hh:64
Data array writers for the VTKWriter.
Definition of the DUNE_DEPRECATED macro for the case that config.h is not available.
A few common exception classes.
Functions for VTK output.
#define DUNE_DEPRECATED_MSG(text)
Mark some entity as deprecated.
Definition: deprecated.hh:169
#define DUNE_THROW(E, m)
Definition: exceptions.hh:243
Traits::Grid Grid
type of the grid
Definition: gridview.hh:77
const CollectiveCommunication & comm() const
obtain collective communication object
Definition: gridview.hh:249
Traits::IndexSet IndexSet
type of the index set
Definition: gridview.hh:80
Grid::ctype ctype
type used for coordinates in grid
Definition: gridview.hh:127
@ dimension
The dimension of the grid.
Definition: gridview.hh:130
@ dimensionworld
The dimension of the world the grid lives in.
Definition: gridview.hh:134
PartitionIteratorType
Parameter to be used for the parallel level- and leaf iterators.
Definition: gridenums.hh:134
PartitionType
Attributes used in the generic overlap model.
Definition: gridenums.hh:28
@ All_Partition
all entities
Definition: gridenums.hh:139
@ InteriorBorder_Partition
interior and border entities
Definition: gridenums.hh:136
@ InteriorEntity
all interior entities
Definition: gridenums.hh:29
std::string relativePath(const std::string &newbase, const std::string &p)
compute a relative path between two paths
Definition: path.cc:151
std::string concatPaths(const std::string &base, const std::string &p)
concatenate two paths
Definition: path.cc:30
Utility class for handling nested indentation in output.
This file implements iterator facade classes for writing stl conformant iterators.
Mapper for multiple codim and multiple geometry types.
Dune namespace.
Definition: alignment.hh:10
STL namespace.
Utilities for handling filesystem paths.
This file implements the class shared_ptr (a reference counting pointer), for those systems that don'...
Static tag representing a codimension.
Definition: dimension.hh:22
static const ReferenceElement< ctype, dim > & general(const GeometryType &type)
get general reference elements
Definition: referenceelements.hh:484
Base class for polymorphic container of underlying data set.
Definition: vtkwriter.hh:157
virtual void write(const Coordinate &pos, Writer &w, std::size_t count) const =0
Evaluate data set at local position pos inside the current entity and write result to w.
virtual void unbind()=0
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
virtual void bind(const Entity &e)=0
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Type erasure implementation for functions conforming to the dune-functions LocalFunction interface.
Definition: vtkwriter.hh:180
virtual void write(const Coordinate &pos, Writer &w, std::size_t count) const
Evaluate data set at local position pos inside the current entity and write result to w.
Definition: vtkwriter.hh:198
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:193
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:188
Type erasure implementation for legacy VTKFunctions.
Definition: vtkwriter.hh:227
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:238
virtual void write(const Coordinate &pos, Writer &w, std::size_t count) const
Evaluate data set at local position pos inside the current entity and write result to w.
Definition: vtkwriter.hh:243
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:233
Definition: typetraits.hh:438
Traits for type conversions and type information.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Dec 22, 23:30, 2024)