dune-grid  2.3.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 
17 #include <dune/common/exceptions.hh>
18 #include <dune/common/indent.hh>
19 #include <dune/common/iteratorfacades.hh>
20 #include <dune/common/path.hh>
21 #include <dune/common/shared_ptr.hh>
22 #include <dune/geometry/referenceelements.hh>
31 
49 namespace Dune
50 {
59  template< class GridView >
60  class VTKWriter {
61  // extract types
62  typedef typename GridView::Grid Grid;
63  typedef typename GridView::ctype DT;
64  enum { n = GridView::dimension };
65  enum { w = GridView::dimensionworld };
66 
67  typedef typename GridView::template Codim< 0 >::Entity Cell;
68  typedef typename GridView::template Codim< n >::Entity Vertex;
69  typedef Cell Entity;
70 
71  typedef typename GridView::IndexSet IndexSet;
72 
73  static const PartitionIteratorType VTK_Partition = InteriorBorder_Partition;
74  //static const PartitionIteratorType VTK_Partition = All_Partition;
75 
76  typedef typename GridView::template Codim< 0 >
77  ::template Partition< VTK_Partition >::Iterator
78  GridCellIterator;
79  typedef typename GridView::template Codim< n >
80  ::template Partition< VTK_Partition >::Iterator
81  GridVertexIterator;
82 
84 
85  // return true if entity should be skipped in Vertex and Corner iterator
86  static bool skipEntity( const PartitionType entityType )
87  {
88  switch( VTK_Partition )
89  {
90  // for All_Partition no entity has to be skipped
91  case All_Partition: return false;
92  case InteriorBorder_Partition: return ( entityType != InteriorEntity );
93  default: DUNE_THROW(NotImplemented,"Add check for this partition type");
94  }
95  return false ;
96  }
97 
98  public:
100  typedef shared_ptr< const VTKFunction > VTKFunctionPtr;
101 
102  protected:
103  typedef typename std::list<VTKFunctionPtr>::const_iterator FunctionIterator;
104 
106 
111  class CellIterator : public GridCellIterator
112  {
113  public:
115  CellIterator(const GridCellIterator & x) : GridCellIterator(x) {}
118  const FieldVector<DT,n> position() const
119  {
120  return ReferenceElements<DT,n>::general((*this)->type()).position(0,0);
121  }
122  };
123 
125  {
126  return gridView_.template begin< 0, VTK_Partition >();
127  }
128 
130  {
131  return gridView_.template end< 0, VTK_Partition >();
132  }
133 
135 
150  public ForwardIteratorFacade<VertexIterator, const Entity, const Entity&, int>
151  {
152  GridCellIterator git;
153  GridCellIterator gend;
154  VTK::DataMode datamode;
155  // Index of the currently visited corner within the current element.
156  // NOTE: this is in Dune-numbering, in contrast to CornerIterator.
157  int cornerIndexDune;
158  const VertexMapper & vertexmapper;
159  std::vector<bool> visited;
160  // in conforming mode, for each vertex id (as obtained by vertexmapper)
161  // hold its number in the iteration order (VertexIterator)
162  int offset;
163  protected:
165  {
166  if( git == gend )
167  return;
168  ++cornerIndexDune;
169  const int numCorners = git->template count< n >();
170  if( cornerIndexDune == numCorners )
171  {
172  offset += numCorners;
173  cornerIndexDune = 0;
174 
175  ++git;
176  while( (git != gend) && skipEntity( git->partitionType() ) )
177  ++git;
178  }
179  }
180  public:
181  VertexIterator(const GridCellIterator & x,
182  const GridCellIterator & end,
183  const VTK::DataMode & dm,
184  const VertexMapper & vm) :
185  git(x), gend(end), datamode(dm), cornerIndexDune(0),
186  vertexmapper(vm), visited(vm.size(), false),
187  offset(0)
188  {
189  if (datamode == VTK::conforming && git != gend)
190  visited[vertexmapper.map(*git,cornerIndexDune,n)] = true;
191  }
192  void increment ()
193  {
194  switch (datamode)
195  {
196  case VTK::conforming :
197  while(visited[vertexmapper.map(*git,cornerIndexDune,n)])
198  {
199  basicIncrement();
200  if (git == gend) return;
201  }
202  visited[vertexmapper.map(*git,cornerIndexDune,n)] = true;
203  break;
204  case VTK::nonconforming :
205  basicIncrement();
206  break;
207  }
208  }
209  bool equals (const VertexIterator & cit) const
210  {
211  return git == cit.git
212  && cornerIndexDune == cit.cornerIndexDune
213  && datamode == cit.datamode;
214  }
215  const Entity& dereference() const
216  {
217  return *git;
218  }
220  int localindex () const
221  {
222  return cornerIndexDune;
223  }
225  const FieldVector<DT,n> & position () const
226  {
227  return ReferenceElements<DT,n>::general(git->type())
228  .position(cornerIndexDune,n);
229  }
230  };
231 
233  {
234  return VertexIterator( gridView_.template begin< 0, VTK_Partition >(),
235  gridView_.template end< 0, VTK_Partition >(),
236  datamode, *vertexmapper );
237  }
238 
240  {
241  return VertexIterator( gridView_.template end< 0, VTK_Partition >(),
242  gridView_.template end< 0, VTK_Partition >(),
243  datamode, *vertexmapper );
244  }
245 
247 
262  public ForwardIteratorFacade<CornerIterator, const Entity, const Entity&, int>
263  {
264  GridCellIterator git;
265  GridCellIterator gend;
266  VTK::DataMode datamode;
267  // Index of the currently visited corner within the current element.
268  // NOTE: this is in VTK-numbering, in contrast to VertexIterator.
269  int cornerIndexVTK;
270  const VertexMapper & vertexmapper;
271  // in conforming mode, for each vertex id (as obtained by vertexmapper)
272  // hold its number in the iteration order of VertexIterator (*not*
273  // CornerIterator)
274  const std::vector<int> & number;
275  // holds the number of corners of all the elements we have seen so far,
276  // excluding the current element
277  int offset;
278 
279  public:
280  CornerIterator(const GridCellIterator & x,
281  const GridCellIterator & end,
282  const VTK::DataMode & dm,
283  const VertexMapper & vm,
284  const std::vector<int> & num) :
285  git(x), gend(end), datamode(dm), cornerIndexVTK(0),
286  vertexmapper(vm),
287  number(num), offset(0) {}
288  void increment ()
289  {
290  if( git == gend )
291  return;
292  ++cornerIndexVTK;
293  const int numCorners = git->template count< n >();
294  if( cornerIndexVTK == numCorners )
295  {
296  offset += numCorners;
297  cornerIndexVTK = 0;
298 
299  ++git;
300  while( (git != gend) && skipEntity( git->partitionType() ) )
301  ++git;
302  }
303  }
304  bool equals (const CornerIterator & cit) const
305  {
306  return git == cit.git
307  && cornerIndexVTK == cit.cornerIndexVTK
308  && datamode == cit.datamode;
309  }
310  const Entity& dereference() const
311  {
312  return *git;
313  }
315 
319  int id () const
320  {
321  switch (datamode)
322  {
323  case VTK::conforming :
324  return
325  number[vertexmapper.map(*git,VTK::renumber(*git,cornerIndexVTK),
326  n)];
327  case VTK::nonconforming :
328  return offset + VTK::renumber(*git,cornerIndexVTK);
329  default :
330  DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
331  }
332  }
333  };
334 
336  {
337  return CornerIterator( gridView_.template begin< 0, VTK_Partition >(),
338  gridView_.template end< 0, VTK_Partition >(),
339  datamode, *vertexmapper, number );
340  }
341 
343  {
344  return CornerIterator( gridView_.template end< 0, VTK_Partition >(),
345  gridView_.template end< 0, VTK_Partition >(),
346  datamode, *vertexmapper, number );
347  }
348 
349  public:
357  explicit VTKWriter ( const GridView &gridView,
359  : gridView_( gridView ),
360  datamode( dm )
361  { }
362 
367  void addCellData (const VTKFunctionPtr & p)
368  {
369  celldata.push_back(p);
370  }
371 
378  {
379  celldata.push_back(VTKFunctionPtr(p));
380  }
381 
397  template<class V>
398  void addCellData (const V& v, const std::string &name, int ncomps = 1)
399  {
400  typedef P0VTKFunction<GridView, V> Function;
401  for (int c=0; c<ncomps; ++c) {
402  std::stringstream compName;
403  compName << name;
404  if (ncomps>1)
405  compName << "[" << c << "]";
406  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
407  celldata.push_back(VTKFunctionPtr(p));
408  }
409  }
410 
417  {
418  vertexdata.push_back(VTKFunctionPtr(p));
419  }
420 
425  void addVertexData (const VTKFunctionPtr & p)
426  {
427  vertexdata.push_back(p);
428  }
429 
445  template<class V>
446  void addVertexData (const V& v, const std::string &name, int ncomps=1)
447  {
448  typedef P1VTKFunction<GridView, V> Function;
449  for (int c=0; c<ncomps; ++c) {
450  std::stringstream compName;
451  compName << name;
452  if (ncomps>1)
453  compName << "[" << c << "]";
454  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
455  vertexdata.push_back(VTKFunctionPtr(p));
456  }
457  }
458 
460  void clear ()
461  {
462  celldata.clear();
463  vertexdata.clear();
464  }
465 
467  virtual ~VTKWriter ()
468  {
469  this->clear();
470  }
471 
483  std::string write ( const std::string &name,
484  VTK::OutputType type = VTK::ascii )
485  {
486  return write( name, type, gridView_.comm().rank(), gridView_.comm().size() );
487  }
488 
515  std::string pwrite ( const std::string & name, const std::string & path, const std::string & extendpath,
516  VTK::OutputType type = VTK::ascii )
517  {
518  return pwrite( name, path, extendpath, type, gridView_.comm().rank(), gridView_.comm().size() );
519  }
520 
521  protected:
523 
534  std::string getParallelPieceName(const std::string& name,
535  const std::string& path,
536  int commRank, int commSize) const
537  {
538  std::ostringstream s;
539  if(path.size() > 0) {
540  s << path;
541  if(path[path.size()-1] != '/')
542  s << '/';
543  }
544  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
545  s << 'p' << std::setw(4) << std::setfill('0') << commRank << '-';
546  s << name;
547  if(GridView::dimension > 1)
548  s << ".vtu";
549  else
550  s << ".vtp";
551  return s.str();
552  }
553 
555 
565  std::string getParallelHeaderName(const std::string& name,
566  const std::string& path,
567  int commSize) const
568  {
569  std::ostringstream s;
570  if(path.size() > 0) {
571  s << path;
572  if(path[path.size()-1] != '/')
573  s << '/';
574  }
575  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
576  s << name;
577  if(GridView::dimension > 1)
578  s << ".pvtu";
579  else
580  s << ".pvtp";
581  return s.str();
582  }
583 
585 
597  std::string getSerialPieceName(const std::string& name,
598  const std::string& path) const
599  {
600  static const std::string extension =
601  GridView::dimension == 1 ? ".vtp" : ".vtu";
602 
603  return concatPaths(path, name+extension);
604  }
605 
621  std::string write ( const std::string &name,
622  VTK::OutputType type,
623  const int commRank,
624  const int commSize )
625  {
626  // in the parallel case, just use pwrite, it has all the necessary
627  // stuff, so we don't need to reimplement it here.
628  if(commSize > 1)
629  return pwrite(name, "", "", type, commRank, commSize);
630 
631  // make data mode visible to private functions
632  outputtype = type;
633 
634  // generate filename for process data
635  std::string pieceName = getSerialPieceName(name, "");
636 
637  // write process data
638  std::ofstream file;
639  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
640  std::ios_base::eofbit);
641  file.open( pieceName.c_str(), std::ios::binary );
642  if (! file.is_open())
643  DUNE_THROW(IOError, "Could not write to piece file " << pieceName);
644  writeDataFile( file );
645  file.close();
646 
647  return pieceName;
648  }
649 
651 
674  std::string pwrite(const std::string& name, const std::string& path,
675  const std::string& extendpath,
676  VTK::OutputType ot, const int commRank,
677  const int commSize )
678  {
679  // make data mode visible to private functions
680  outputtype=ot;
681 
682  // do some magic because paraview can only cope with relative pathes to piece files
683  std::ofstream file;
684  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
685  std::ios_base::eofbit);
686  std::string piecepath = concatPaths(path, extendpath);
687  std::string relpiecepath = relativePath(path, piecepath);
688 
689  // write this processes .vtu/.vtp piece file
690  std::string fullname = getParallelPieceName(name, piecepath, commRank,
691  commSize);
692  file.open(fullname.c_str(),std::ios::binary);
693  if (! file.is_open())
694  DUNE_THROW(IOError, "Could not write to piecefile file " << fullname);
695  writeDataFile(file);
696  file.close();
697  gridView_.comm().barrier();
698 
699  // if we are rank 0, write .pvtu/.pvtp parallel header
700  fullname = getParallelHeaderName(name, path, commSize);
701  if( commRank ==0 )
702  {
703  file.open(fullname.c_str());
704  if (! file.is_open())
705  DUNE_THROW(IOError, "Could not write to parallel file " << fullname);
706  writeParallelHeader(file,name,relpiecepath, commSize );
707  file.close();
708  }
709  gridView_.comm().barrier();
710  return fullname;
711  }
712 
713  private:
715 
732  void writeParallelHeader(std::ostream& s, const std::string& piecename,
733  const std::string& piecepath, const int commSize)
734  {
735  VTK::FileType fileType =
737 
738  VTK::PVTUWriter writer(s, fileType);
739 
740  writer.beginMain();
741 
742  // PPointData
743  {
744  std::string scalars;
745  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end();
746  ++it)
747  if ((*it)->ncomps()==1)
748  {
749  scalars = (*it)->name();
750  break;
751  }
752  std::string vectors;
753  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end();
754  ++it)
755  if ((*it)->ncomps()>1)
756  {
757  vectors = (*it)->name();
758  break;
759  }
760  writer.beginPointData(scalars, vectors);
761  }
762  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end();
763  ++it)
764  {
765  unsigned writecomps = (*it)->ncomps();
766  if(writecomps == 2) writecomps = 3;
767  writer.addArray<float>((*it)->name(), writecomps);
768  }
769  writer.endPointData();
770 
771  // PCellData
772  {
773  std::string scalars;
774  for (FunctionIterator it=celldata.begin(); it!=celldata.end();
775  ++it)
776  if ((*it)->ncomps()==1)
777  {
778  scalars = (*it)->name();
779  break;
780  }
781  std::string vectors;
782  for (FunctionIterator it=celldata.begin(); it!=celldata.end();
783  ++it)
784  if ((*it)->ncomps()>1)
785  {
786  vectors = (*it)->name();
787  break;
788  }
789  writer.beginCellData(scalars, vectors);
790  }
791  for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it) {
792  unsigned writecomps = (*it)->ncomps();
793  if(writecomps == 2) writecomps = 3;
794  writer.addArray<float>((*it)->name(), writecomps);
795  }
796  writer.endCellData();
797 
798  // PPoints
799  writer.beginPoints();
800  writer.addArray<float>("Coordinates", 3);
801  writer.endPoints();
802 
803  // Pieces
804  for( int i = 0; i < commSize; ++i )
805  {
806  const std::string& fullname = getParallelPieceName(piecename,
807  piecepath, i,
808  commSize);
809  writer.addPiece(fullname);
810  }
811 
812  writer.endMain();
813  }
814 
816  void writeDataFile (std::ostream& s)
817  {
818  VTK::FileType fileType =
820 
821  VTK::VTUWriter writer(s, outputtype, fileType);
822 
823  // Grid characteristics
824  vertexmapper = new VertexMapper( gridView_ );
825  if (datamode == VTK::conforming)
826  {
827  number.resize(vertexmapper->size());
828  for (std::vector<int>::size_type i=0; i<number.size(); i++) number[i] = -1;
829  }
831 
832  writer.beginMain(ncells, nvertices);
833  writeAllData(writer);
834  writer.endMain();
835 
836  // write appended binary data section
837  if(writer.beginAppended())
838  writeAllData(writer);
839  writer.endAppended();
840 
841  delete vertexmapper; number.clear();
842  }
843 
844  void writeAllData(VTK::VTUWriter& writer) {
845  // PointData
846  writeVertexData(writer);
847 
848  // CellData
849  writeCellData(writer);
850 
851  // Points
852  writeGridPoints(writer);
853 
854  // Cells
855  writeGridCells(writer);
856  }
857 
858  protected:
859  std::string getFormatString() const
860  {
861  if (outputtype==VTK::ascii)
862  return "ascii";
863  if (outputtype==VTK::base64)
864  return "binary";
866  return "appended";
868  return "appended";
869  DUNE_THROW(IOError, "VTKWriter: unsupported OutputType" << outputtype);
870  }
871 
872  std::string getTypeString() const
873  {
874  if (n==1)
875  return "PolyData";
876  else
877  return "UnstructuredGrid";
878  }
879 
881  virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
882  {
883  nvertices = 0;
884  ncells = 0;
885  ncorners = 0;
886  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
887  {
888  ncells++;
889  // because of the use of vertexmapper->map(), this iteration must be
890  // in the order of Dune's numbering.
891  for (int i=0; i<it->template count<n>(); ++i)
892  {
893  ncorners++;
894  if (datamode == VTK::conforming)
895  {
896  int alpha = vertexmapper->map(*it,i,n);
897  if (number[alpha]<0)
898  number[alpha] = nvertices++;
899  }
900  else
901  {
902  nvertices++;
903  }
904  }
905  }
906  }
907 
909  virtual void writeCellData(VTK::VTUWriter& writer)
910  {
911  if(celldata.size() == 0)
912  return;
913 
914  std::string scalars = "";
915  for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
916  if ((*it)->ncomps()==1)
917  {
918  scalars = (*it)->name();
919  break;
920  }
921  std::string vectors = "";
922  for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
923  if ((*it)->ncomps()>1)
924  {
925  vectors = (*it)->name();
926  break;
927  }
928 
929  writer.beginCellData(scalars, vectors);
930  for (FunctionIterator it=celldata.begin(); it!=celldata.end(); ++it)
931  {
932  // vtk file format: a vector data always should have 3 comps (with
933  // 3rd comp = 0 in 2D case)
934  unsigned writecomps = (*it)->ncomps();
935  if(writecomps == 2) writecomps = 3;
936  shared_ptr<VTK::DataArrayWriter<float> > p
937  (writer.makeArrayWriter<float>((*it)->name(), writecomps,
938  ncells));
939  if(!p->writeIsNoop())
940  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
941  {
942  for (int j=0; j<(*it)->ncomps(); j++)
943  p->write((*it)->evaluate(j,*i,i.position()));
944  // vtk file format: a vector data always should have 3 comps
945  // (with 3rd comp = 0 in 2D case)
946  for (unsigned j=(*it)->ncomps(); j < writecomps; ++j)
947  p->write(0.0);
948  }
949  }
950  writer.endCellData();
951  }
952 
954  virtual void writeVertexData(VTK::VTUWriter& writer)
955  {
956  if(vertexdata.size() == 0)
957  return;
958 
959  std::string scalars = "";
960  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
961  if ((*it)->ncomps()==1)
962  {
963  scalars = (*it)->name();
964  break;
965  }
966  std::string vectors = "";
967  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
968  if ((*it)->ncomps()>1)
969  {
970  vectors = (*it)->name();
971  break;
972  }
973 
974  writer.beginPointData(scalars, vectors);
975  for (FunctionIterator it=vertexdata.begin(); it!=vertexdata.end(); ++it)
976  {
977  // vtk file format: a vector data always should have 3 comps (with
978  // 3rd comp = 0 in 2D case)
979  unsigned writecomps = (*it)->ncomps();
980  if(writecomps == 2) writecomps = 3;
981  shared_ptr<VTK::DataArrayWriter<float> > p
982  (writer.makeArrayWriter<float>((*it)->name(), writecomps,
983  nvertices));
984  if(!p->writeIsNoop())
985  for (VertexIterator vit=vertexBegin(); vit!=vertexEnd(); ++vit)
986  {
987  for (int j=0; j<(*it)->ncomps(); j++)
988  p->write((*it)->evaluate(j,*vit,vit.position()));
989  // vtk file format: a vector data always should have 3 comps
990  // (with 3rd comp = 0 in 2D case)
991  for (unsigned j=(*it)->ncomps(); j < writecomps; ++j)
992  p->write(0.0);
993  }
994  }
995  writer.endPointData();
996  }
997 
999  virtual void writeGridPoints(VTK::VTUWriter& writer)
1000  {
1001  writer.beginPoints();
1002 
1003  shared_ptr<VTK::DataArrayWriter<float> > p
1004  (writer.makeArrayWriter<float>("Coordinates", 3, nvertices));
1005  if(!p->writeIsNoop()) {
1006  VertexIterator vEnd = vertexEnd();
1007  for (VertexIterator vit=vertexBegin(); vit!=vEnd; ++vit)
1008  {
1009  int dimw=w;
1010  for (int j=0; j<std::min(dimw,3); j++)
1011  p->write(vit->geometry().corner(vit.localindex())[j]);
1012  for (int j=std::min(dimw,3); j<3; j++)
1013  p->write(0.0);
1014  }
1015  }
1016  // free the VTK::DataArrayWriter before touching the stream
1017  p.reset();
1018 
1019  writer.endPoints();
1020  }
1021 
1023  virtual void writeGridCells(VTK::VTUWriter& writer)
1024  {
1025  writer.beginCells();
1026 
1027  // connectivity
1028  {
1029  shared_ptr<VTK::DataArrayWriter<int> > p1
1030  (writer.makeArrayWriter<int>("connectivity", 1, ncorners));
1031  if(!p1->writeIsNoop())
1032  for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
1033  p1->write(it.id());
1034  }
1035 
1036  // offsets
1037  {
1038  shared_ptr<VTK::DataArrayWriter<int> > p2
1039  (writer.makeArrayWriter<int>("offsets", 1, ncells));
1040  if(!p2->writeIsNoop()) {
1041  int offset = 0;
1042  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1043  {
1044  offset += it->template count<n>();
1045  p2->write(offset);
1046  }
1047  }
1048  }
1049 
1050  // types
1051  if (n>1)
1052  {
1053  shared_ptr<VTK::DataArrayWriter<unsigned char> > p3
1054  (writer.makeArrayWriter<unsigned char>("types", 1, ncells));
1055  if(!p3->writeIsNoop())
1056  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1057  {
1058  int vtktype = VTK::geometryType(it->type());
1059  p3->write(vtktype);
1060  }
1061  }
1062 
1063  writer.endCells();
1064  }
1065 
1066  protected:
1067  // the list of registered functions
1068  std::list<VTKFunctionPtr> celldata;
1069  std::list<VTKFunctionPtr> vertexdata;
1070 
1071  // the grid
1073 
1074  // temporary grid information
1075  int ncells;
1078  private:
1079  VertexMapper* vertexmapper;
1080  // in conforming mode, for each vertex id (as obtained by vertexmapper)
1081  // hold its number in the iteration order (VertexIterator)
1082  std::vector<int> number;
1083  VTK::DataMode datamode;
1084  protected:
1086  };
1087 
1088 }
1089 
1090 #endif
Traits::IndexSet IndexSet
type of the index set
Definition: common/gridview.hh:70
const Entity & dereference() const
Definition: vtkwriter.hh:310
The dimension of the world the grid lives in.
Definition: common/gridview.hh:124
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: vtkwriter.hh:999
GridView gridView_
Definition: vtkwriter.hh:1072
PartitionIteratorType
Parameter to be used for the parallel level- and leaf iterators.
Definition: gridenums.hh:130
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:621
std::string getFormatString() const
Definition: vtkwriter.hh:859
void increment()
Definition: vtkwriter.hh:288
CornerIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm, const std::vector< int > &num)
Definition: vtkwriter.hh:280
CornerIterator cornerBegin() const
Definition: vtkwriter.hh:335
Data array writers for the VTKWriter.
Writer for the ouput of grid functions in the vtk format.Writes arbitrary grid functions (living on c...
Definition: vtkwriter.hh:60
Output to the file is inline base64 binary.
Definition: common.hh:44
FileType
which type of VTK file to write
Definition: common.hh:290
Traits::Grid Grid
type of the grid
Definition: common/gridview.hh:67
VertexIterator vertexBegin() const
Definition: vtkwriter.hh:232
Implementation class for a multiple codim and multiple geometry type mapper.
Definition: mcmgmapper.hh:102
std::list< VTKFunctionPtr > celldata
Definition: vtkwriter.hh:1068
all entities
Definition: gridenums.hh:135
for .vtu files (UnstructuredGrid)
Definition: common.hh:294
The dimension of the grid.
Definition: common/gridview.hh:120
for .vtp files (PolyData)
Definition: common.hh:292
Grid::ctype ctype
type used for coordinates in grid
Definition: common/gridview.hh:117
void endCellData()
finish CellData section
Definition: vtuwriter.hh:218
A base class for grid functions with any return type and dimension.
Definition: function.hh:37
Common stuff for the VTKWriter.
void basicIncrement()
Definition: vtkwriter.hh:164
Take a vector and interpret it as cell data for the VTKWriter.
Definition: function.hh:87
Output to the file is in ascii.
Definition: common.hh:42
void increment()
Definition: vtkwriter.hh:192
void beginPoints()
start section for the point coordinates
Definition: vtuwriter.hh:236
int size() const
Return total number of entities in the entity set managed by the mapper.
Definition: mcmgmapper.hh:181
Functions for VTK output.
Iterate over the elements' corners.
Definition: vtkwriter.hh:261
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: vtkwriter.hh:909
const Entity & dereference() const
Definition: vtkwriter.hh:215
void clear()
clear list of registered functions
Definition: vtkwriter.hh:460
Iterate over the grid's vertices.
Definition: vtkwriter.hh:149
Mapper for multiple codim and multiple geometry types.
int localindex() const
index of vertex within the entity, in Dune-numbering
Definition: vtkwriter.hh:220
std::string getTypeString() const
Definition: vtkwriter.hh:872
CornerIterator cornerEnd() const
Definition: vtkwriter.hh:342
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: vtkwriter.hh:954
int ncells
Definition: vtkwriter.hh:1075
void addVertexData(const VTKFunctionPtr &p)
Add a grid function that lives on the vertices of the grid to the visualization.
Definition: vtkwriter.hh:425
CellIterator cellEnd() const
Definition: vtkwriter.hh:129
int ncorners
Definition: vtkwriter.hh:1077
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: vtkwriter.hh:1023
GeometryType geometryType(const Dune::GeometryType &t)
mapping from GeometryType to VTKGeometryType
Definition: common.hh:195
Iterator over the grids elements.
Definition: vtkwriter.hh:111
std::list< VTKFunctionPtr > vertexdata
Definition: vtkwriter.hh:1069
Output non-conforming data.
Definition: common.hh:78
void endPointData()
finish PointData section
Definition: vtuwriter.hh:180
Grid view abstract base classInterface class for a view on grids. Grids return two types of view...
Definition: common/gridview.hh:56
CellIterator(const GridCellIterator &x)
construct a CellIterator from the gridview's Iterator.
Definition: vtkwriter.hh:115
PartitionType
Attributes used in the generic overlap model.
Definition: gridenums.hh:24
shared_ptr< const VTKFunction > VTKFunctionPtr
Definition: vtkwriter.hh:100
Dump a .vtu/.vtp files contents to a stream.
Definition: pvtuwriter.hh:60
void addCellData(const VTKFunctionPtr &p)
Add a grid function that lives on the cells of the grid to the visualization.
Definition: vtkwriter.hh:367
CellIterator cellBegin() const
Definition: vtkwriter.hh:124
Ouput is to the file is appended raw binary.
Definition: common.hh:46
const FieldVector< DT, n > & position() const
position of vertex inside the entity
Definition: vtkwriter.hh:225
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:40
void endCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:283
VTK::OutputType outputtype
Definition: vtkwriter.hh:1085
std::string getSerialPieceName(const std::string &name, const std::string &path) const
return name of a serial piece file
Definition: vtkwriter.hh:597
const FieldVector< DT, n > position() const
Definition: vtkwriter.hh:118
bool equals(const CornerIterator &cit) const
Definition: vtkwriter.hh:304
VertexIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm)
Definition: vtkwriter.hh:181
Output conforming data.
Definition: common.hh:70
interior and border entities
Definition: gridenums.hh:132
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:398
void beginCellData(const std::string &scalars="", const std::string &vectors="")
start CellData section
Definition: vtuwriter.hh:203
int id() const
Process-local consecutive zero-starting vertex id.
Definition: vtkwriter.hh:319
VTKWriter(const GridView &gridView, VTK::DataMode dm=VTK::conforming)
Construct a VTKWriter working on a specific GridView.
Definition: vtkwriter.hh:357
virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
count the vertices, cells and corners
Definition: vtkwriter.hh:881
Dune::VTKFunction< GridView > VTKFunction
Definition: vtkwriter.hh:99
DataMode
Whether to produce conforming or non-conforming output.
Definition: common.hh:64
int renumber(const Dune::GeometryType &t, int i)
renumber VTK <-> Dune
Definition: common.hh:224
int map(const EntityType &e) const
Map entity to array index.
Definition: mcmgmapper.hh:153
Ouput is to the file is appended base64 binary.
Definition: common.hh:48
int nvertices
Definition: vtkwriter.hh:1076
std::string write(const std::string &name, VTK::OutputType type=VTK::ascii)
write output (interface might change later)
Definition: vtkwriter.hh:483
bool equals(const VertexIterator &cit) const
Definition: vtkwriter.hh:209
void addCellData(VTKFunction *p)
Add a grid function that lives on the cells of the grid to the visualization.
Definition: vtkwriter.hh:377
std::string getParallelHeaderName(const std::string &name, const std::string &path, int commSize) const
return name of a parallel header file
Definition: vtkwriter.hh:565
const CollectiveCommunication & comm() const
obtain collective communication object
Definition: common/gridview.hh:232
void beginCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:272
void addVertexData(VTKFunction *p)
Add a grid function that lives on the vertices of the grid to the visualization.
Definition: vtkwriter.hh:416
Dump a .vtu/.vtp files contents to a stream.
Definition: vtuwriter.hh:96
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:515
void endPoints()
finish section for the point coordinates
Definition: vtuwriter.hh:247
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:534
void beginPointData(const std::string &scalars="", const std::string &vectors="")
start PointData section
Definition: vtuwriter.hh:165
VertexIterator vertexEnd() const
Definition: vtkwriter.hh:239
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:674
virtual ~VTKWriter()
destructor
Definition: vtkwriter.hh:467
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:446
all interior entities
Definition: gridenums.hh:25
int min(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:352
DataArrayWriter< T > * makeArrayWriter(const std::string &name, unsigned ncomps, unsigned nitems)
aquire a DataArrayWriter
Definition: vtuwriter.hh:379
Take a vector and interpret it as point data for the VTKWriter.
Definition: function.hh:190
std::list< VTKFunctionPtr >::const_iterator FunctionIterator
Definition: vtkwriter.hh:103