dune-grid  3.0-git
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 #include <memory>
14 
15 #include <vector>
16 #include <list>
17 
18 #include <dune/common/typetraits.hh>
19 #include <dune/common/exceptions.hh>
20 #include <dune/common/std/memory.hh>
21 #include <dune/common/indent.hh>
22 #include <dune/common/iteratorfacades.hh>
23 #include <dune/common/path.hh>
24 #include <dune/geometry/referenceelements.hh>
33 
47 namespace Dune
48 {
49 
50  namespace detail {
51 
52  template<typename F, typename = int>
54  : public std::false_type
55  {};
56 
57  template<typename T>
58  struct _has_local_context<T,typename std::enable_if<(sizeof(std::declval<T>().localContext()) > 0),int>::type>
59  : public std::true_type
60  {};
61 
62  }
63 
64  namespace VTKWriteTypeTraits {
65  template<typename T>
67  {
68  };
69  }
70 
71  // Forward-declaration here, so the class can be friend of VTKWriter
72  template <class GridView>
74  template <class GridView>
75  class VTKSequenceWriter;
76 
85  template< class GridView >
86  class VTKWriter {
87 
88  // VTKSequenceWriterBase needs getSerialPieceName
89  // and getParallelHeaderName
91  // VTKSequenceWriter needs the grid view, to get the MPI size and rank
92  friend class VTKSequenceWriter<GridView>;
93 
94  // extract types
95  typedef typename GridView::Grid Grid;
96  typedef typename GridView::ctype DT;
97  enum { n = GridView::dimension };
98  enum { w = GridView::dimensionworld };
99 
100  typedef typename GridView::template Codim< 0 >::Entity Cell;
101  typedef typename GridView::template Codim< n >::Entity Vertex;
102  typedef Cell Entity;
103 
104  typedef typename GridView::IndexSet IndexSet;
105 
106  static const PartitionIteratorType VTK_Partition = InteriorBorder_Partition;
107  //static const PartitionIteratorType VTK_Partition = All_Partition;
108 
109  typedef typename GridView::template Codim< 0 >
110  ::template Partition< VTK_Partition >::Iterator
111  GridCellIterator;
112  typedef typename GridView::template Codim< n >
113  ::template Partition< VTK_Partition >::Iterator
114  GridVertexIterator;
115 
116  typedef typename GridCellIterator::Reference EntityReference;
117 
118  typedef typename GridView::template Codim< 0 >
119  ::Entity::Geometry::LocalCoordinate Coordinate;
120 
121  typedef MultipleCodimMultipleGeomTypeMapper< GridView, MCMGVertexLayout > VertexMapper;
122 
123  // return true if entity should be skipped in Vertex and Corner iterator
124  static bool skipEntity( const PartitionType entityType )
125  {
126  switch( VTK_Partition )
127  {
128  // for All_Partition no entity has to be skipped
129  case All_Partition: return false;
130  case InteriorBorder_Partition: return ( entityType != InteriorEntity );
131  default: DUNE_THROW(NotImplemented,"Add check for this partition type");
132  }
133  return false ;
134  }
135 
136  public:
137 
139 
140  protected:
141 
143 
147  {
148 
149  public:
150 
152 
155  {
156 
158  virtual void bind(const Entity& e) = 0;
159 
161  virtual void unbind() = 0;
162 
164 
167  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const = 0;
168 
170  {}
171 
172  };
173 
175  template<typename F>
177  : public FunctionWrapperBase
178  {
179  using Function = typename std::decay<F>::type;
180 
181  template<typename F_>
183  : _f(std::forward<F_>(f))
184  {}
185 
186  virtual void bind(const Entity& e)
187  {
188  _f.bind(e);
189  }
190 
191  virtual void unbind()
192  {
193  _f.unbind();
194  }
195 
196  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
197  {
198  auto r = _f(pos);
199  // we need to do different things here depending on whether r supports indexing into it or not.
200  do_write(w,r,count,is_indexable<decltype(r)>());
201  }
202 
203  private:
204 
205  template<typename R>
206  void do_write(Writer& w, const R& r, std::size_t count, std::true_type) const
207  {
208  for (std::size_t i = 0; i < count; ++i)
209  w.write(r[i]);
210  }
211 
212  template<typename R>
213  void do_write(Writer& w, const R& r, std::size_t count, std::false_type) const
214  {
215  assert(count == 1);
216  w.write(r);
217  }
218 
219  Function _f;
220  };
221 
224  : public FunctionWrapperBase
225  {
226  VTKFunctionWrapper(const std::shared_ptr< const VTKFunction >& f)
227  : _f(f)
228  , _entity(nullptr)
229  {}
230 
231  virtual void bind(const Entity& e)
232  {
233  _entity = &e;
234  }
235 
236  virtual void unbind()
237  {
238  _entity = nullptr;
239  }
240 
241  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
242  {
243  for (std::size_t i = 0; i < count; ++i)
244  w.write(_f->evaluate(i,*_entity,pos));
245  }
246 
247  private:
248 
249  std::shared_ptr< const VTKFunction > _f;
250  const Entity* _entity;
251 
252  };
253 
255  template<typename F>
257  typename std::enable_if<detail::_has_local_context<F>::value,int>::type dummy = 0)
258  : _f(Dune::Std::make_unique<FunctionWrapper<F> >(std::forward<F>(f)))
260  {}
261 
263  template<typename F>
265  typename std::enable_if<not detail::_has_local_context<F>::value,int>::type dummy = 0)
266  : _f(Dune::Std::make_unique< FunctionWrapper<
267  typename std::decay<decltype(localFunction(std::forward<F>(f)))>::type
268  > >(localFunction(std::forward<F>(f))))
270  {}
271 
273  explicit VTKLocalFunction (const std::shared_ptr< const VTKFunction >& vtkFunctionPtr)
274  : _f(Dune::Std::make_unique<VTKFunctionWrapper>(vtkFunctionPtr))
275  , _fieldInfo(
276  vtkFunctionPtr->name(),
277  vtkFunctionPtr->ncomps() > 1 ? VTK::FieldInfo::Type::vector : VTK::FieldInfo::Type::scalar,
278  vtkFunctionPtr->ncomps()
279  )
280  {}
281 
283  std::string name() const
284  {
285  return fieldInfo().name();
286  }
287 
289  const VTK::FieldInfo& fieldInfo() const
290  {
291  return _fieldInfo;
292  }
293 
295  void bind(const Entity& e) const
296  {
297  _f->bind(e);
298  }
299 
301  void unbind() const
302  {
303  _f->unbind();
304  }
305 
307  void write(const Coordinate& pos, Writer& w) const
308  {
309  _f->write(pos,w,fieldInfo().size());
310  }
311 
312  std::shared_ptr<FunctionWrapperBase> _f;
314 
315  };
316 
317  typedef typename std::list<VTKLocalFunction>::const_iterator FunctionIterator;
318 
320 
325  class CellIterator : public GridCellIterator
326  {
327  public:
329  CellIterator(const GridCellIterator & x) : GridCellIterator(x) {}
332  const FieldVector<DT,n> position() const
333  {
334  return ReferenceElements<DT,n>::general((*this)->type()).position(0,0);
335  }
336  };
337 
338  CellIterator cellBegin() const
339  {
340  return gridView_.template begin< 0, VTK_Partition >();
341  }
342 
343  CellIterator cellEnd() const
344  {
345  return gridView_.template end< 0, VTK_Partition >();
346  }
347 
349 
364  public ForwardIteratorFacade<VertexIterator, const Entity, EntityReference, int>
365  {
366  GridCellIterator git;
367  GridCellIterator gend;
368  VTK::DataMode datamode;
369  // Index of the currently visited corner within the current element.
370  // NOTE: this is in Dune-numbering, in contrast to CornerIterator.
371  int cornerIndexDune;
372  const VertexMapper & vertexmapper;
373  std::vector<bool> visited;
374  // in conforming mode, for each vertex id (as obtained by vertexmapper)
375  // hold its number in the iteration order (VertexIterator)
376  int offset;
377 
378  // hide operator ->
379  void operator->();
380  protected:
382  {
383  if( git == gend )
384  return;
385  ++cornerIndexDune;
386  const int numCorners = git->subEntities(n);
387  if( cornerIndexDune == numCorners )
388  {
389  offset += numCorners;
390  cornerIndexDune = 0;
391 
392  ++git;
393  while( (git != gend) && skipEntity( git->partitionType() ) )
394  ++git;
395  }
396  }
397  public:
398  VertexIterator(const GridCellIterator & x,
399  const GridCellIterator & end,
400  const VTK::DataMode & dm,
401  const VertexMapper & vm) :
402  git(x), gend(end), datamode(dm), cornerIndexDune(0),
403  vertexmapper(vm), visited(vm.size(), false),
404  offset(0)
405  {
406  if (datamode == VTK::conforming && git != gend)
407  visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
408  }
409  void increment ()
410  {
411  switch (datamode)
412  {
413  case VTK::conforming :
414  while(visited[vertexmapper.subIndex(*git,cornerIndexDune,n)])
415  {
416  basicIncrement();
417  if (git == gend) return;
418  }
419  visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
420  break;
421  case VTK::nonconforming :
422  basicIncrement();
423  break;
424  }
425  }
426  bool equals (const VertexIterator & cit) const
427  {
428  return git == cit.git
429  && cornerIndexDune == cit.cornerIndexDune
430  && datamode == cit.datamode;
431  }
432  EntityReference dereference() const
433  {
434  return *git;
435  }
437  int localindex () const
438  {
439  return cornerIndexDune;
440  }
442  const FieldVector<DT,n> & position () const
443  {
444  return ReferenceElements<DT,n>::general(git->type())
445  .position(cornerIndexDune,n);
446  }
447  };
448 
449  VertexIterator vertexBegin () const
450  {
451  return VertexIterator( gridView_.template begin< 0, VTK_Partition >(),
452  gridView_.template end< 0, VTK_Partition >(),
453  datamode, *vertexmapper );
454  }
455 
456  VertexIterator vertexEnd () const
457  {
458  return VertexIterator( gridView_.template end< 0, VTK_Partition >(),
459  gridView_.template end< 0, VTK_Partition >(),
460  datamode, *vertexmapper );
461  }
462 
464 
479  public ForwardIteratorFacade<CornerIterator, const Entity, EntityReference, int>
480  {
481  GridCellIterator git;
482  GridCellIterator gend;
483  VTK::DataMode datamode;
484  // Index of the currently visited corner within the current element.
485  // NOTE: this is in VTK-numbering, in contrast to VertexIterator.
486  int cornerIndexVTK;
487  const VertexMapper & vertexmapper;
488  // in conforming mode, for each vertex id (as obtained by vertexmapper)
489  // hold its number in the iteration order of VertexIterator (*not*
490  // CornerIterator)
491  const std::vector<int> & number;
492  // holds the number of corners of all the elements we have seen so far,
493  // excluding the current element
494  int offset;
495 
496  // hide operator ->
497  void operator->();
498  public:
499  CornerIterator(const GridCellIterator & x,
500  const GridCellIterator & end,
501  const VTK::DataMode & dm,
502  const VertexMapper & vm,
503  const std::vector<int> & num) :
504  git(x), gend(end), datamode(dm), cornerIndexVTK(0),
505  vertexmapper(vm),
506  number(num), offset(0) {}
507  void increment ()
508  {
509  if( git == gend )
510  return;
511  ++cornerIndexVTK;
512  const int numCorners = git->subEntities(n);
513  if( cornerIndexVTK == numCorners )
514  {
515  offset += numCorners;
516  cornerIndexVTK = 0;
517 
518  ++git;
519  while( (git != gend) && skipEntity( git->partitionType() ) )
520  ++git;
521  }
522  }
523  bool equals (const CornerIterator & cit) const
524  {
525  return git == cit.git
526  && cornerIndexVTK == cit.cornerIndexVTK
527  && datamode == cit.datamode;
528  }
529  EntityReference dereference() const
530  {
531  return *git;
532  }
534 
538  int id () const
539  {
540  switch (datamode)
541  {
542  case VTK::conforming :
543  return
544  number[vertexmapper.subIndex(*git,VTK::renumber(*git,cornerIndexVTK),
545  n)];
546  case VTK::nonconforming :
547  return offset + VTK::renumber(*git,cornerIndexVTK);
548  default :
549  DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
550  }
551  }
552  };
553 
554  CornerIterator cornerBegin () const
555  {
556  return CornerIterator( gridView_.template begin< 0, VTK_Partition >(),
557  gridView_.template end< 0, VTK_Partition >(),
558  datamode, *vertexmapper, number );
559  }
560 
561  CornerIterator cornerEnd () const
562  {
563  return CornerIterator( gridView_.template end< 0, VTK_Partition >(),
564  gridView_.template end< 0, VTK_Partition >(),
565  datamode, *vertexmapper, number );
566  }
567 
568  public:
576  explicit VTKWriter ( const GridView &gridView,
578  : gridView_( gridView ),
579  datamode( dm )
580  { }
581 
586  void addCellData (const std::shared_ptr< const VTKFunction > & p)
587  {
588  celldata.push_back(VTKLocalFunction(p));
589  }
590 
591  template<typename F>
592  void addCellData(F&& f, VTK::FieldInfo vtkFieldInfo)
593  {
594  celldata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
595  }
596 
612  template<class Container>
613  void addCellData (const Container& v, const std::string &name, int ncomps = 1)
614  {
615  typedef P0VTKFunction<GridView, Container> Function;
616  for (int c=0; c<ncomps; ++c) {
617  std::stringstream compName;
618  compName << name;
619  if (ncomps>1)
620  compName << "[" << c << "]";
621  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
622  addCellData(std::shared_ptr< const VTKFunction >(p));
623  }
624  }
625 
630  void addVertexData (const std::shared_ptr< const VTKFunction > & p)
631  {
632  vertexdata.push_back(VTKLocalFunction(p));
633  }
634 
635  template<typename F>
636  void addVertexData(F&& f, VTK::FieldInfo vtkFieldInfo)
637  {
638  vertexdata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
639  }
640 
641 
657  template<class Container>
658  void addVertexData (const Container& v, const std::string &name, int ncomps=1)
659  {
660  typedef P1VTKFunction<GridView, Container> Function;
661  for (int c=0; c<ncomps; ++c) {
662  std::stringstream compName;
663  compName << name;
664  if (ncomps>1)
665  compName << "[" << c << "]";
666  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c);
667  addVertexData(std::shared_ptr< const VTKFunction >(p));
668  }
669  }
670 
672  void clear ()
673  {
674  celldata.clear();
675  vertexdata.clear();
676  }
677 
679  virtual ~VTKWriter ()
680  {
681  this->clear();
682  }
683 
695  std::string write ( const std::string &name,
696  VTK::OutputType type = VTK::ascii )
697  {
698  return write( name, type, gridView_.comm().rank(), gridView_.comm().size() );
699  }
700 
727  std::string pwrite ( const std::string & name, const std::string & path, const std::string & extendpath,
728  VTK::OutputType type = VTK::ascii )
729  {
730  return pwrite( name, path, extendpath, type, gridView_.comm().rank(), gridView_.comm().size() );
731  }
732 
733  protected:
735 
746  std::string getParallelPieceName(const std::string& name,
747  const std::string& path,
748  int commRank, int commSize) const
749  {
750  std::ostringstream s;
751  if(path.size() > 0) {
752  s << path;
753  if(path[path.size()-1] != '/')
754  s << '/';
755  }
756  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
757  s << 'p' << std::setw(4) << std::setfill('0') << commRank << '-';
758  s << name;
759  if(GridView::dimension > 1)
760  s << ".vtu";
761  else
762  s << ".vtp";
763  return s.str();
764  }
765 
767 
777  std::string getParallelHeaderName(const std::string& name,
778  const std::string& path,
779  int commSize) const
780  {
781  std::ostringstream s;
782  if(path.size() > 0) {
783  s << path;
784  if(path[path.size()-1] != '/')
785  s << '/';
786  }
787  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
788  s << name;
789  if(GridView::dimension > 1)
790  s << ".pvtu";
791  else
792  s << ".pvtp";
793  return s.str();
794  }
795 
797 
809  std::string getSerialPieceName(const std::string& name,
810  const std::string& path) const
811  {
812  static const std::string extension =
813  GridView::dimension == 1 ? ".vtp" : ".vtu";
814 
815  return concatPaths(path, name+extension);
816  }
817 
833  std::string write ( const std::string &name,
834  VTK::OutputType type,
835  const int commRank,
836  const int commSize )
837  {
838  // in the parallel case, just use pwrite, it has all the necessary
839  // stuff, so we don't need to reimplement it here.
840  if(commSize > 1)
841  return pwrite(name, "", "", type, commRank, commSize);
842 
843  // make data mode visible to private functions
844  outputtype = type;
845 
846  // generate filename for process data
847  std::string pieceName = getSerialPieceName(name, "");
848 
849  // write process data
850  std::ofstream file;
851  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
852  std::ios_base::eofbit);
853  // check if file can be opened
854  try {
855  file.open( pieceName.c_str(), std::ios::binary );
856  }
857  catch(...) {
858  std::cerr << "Filename: " << pieceName << " could not be opened" << std::endl;
859  throw;
860  }
861  if (! file.is_open())
862  DUNE_THROW(IOError, "Could not write to piece file " << pieceName);
863  writeDataFile( file );
864  file.close();
865 
866  return pieceName;
867  }
868 
870 
893  std::string pwrite(const std::string& name, const std::string& path,
894  const std::string& extendpath,
895  VTK::OutputType ot, const int commRank,
896  const int commSize )
897  {
898  // make data mode visible to private functions
899  outputtype=ot;
900 
901  // do some magic because paraview can only cope with relative paths to piece files
902  std::ofstream file;
903  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
904  std::ios_base::eofbit);
905  std::string piecepath = concatPaths(path, extendpath);
906  std::string relpiecepath = relativePath(path, piecepath);
907 
908  // write this processes .vtu/.vtp piece file
909  std::string fullname = getParallelPieceName(name, piecepath, commRank,
910  commSize);
911  // check if file can be opened
912  try {
913  file.open(fullname.c_str(),std::ios::binary);
914  }
915  catch(...) {
916  std::cerr << "Filename: " << fullname << " could not be opened" << std::endl;
917  throw;
918  }
919  if (! file.is_open())
920  DUNE_THROW(IOError, "Could not write to piecefile file " << fullname);
921  writeDataFile(file);
922  file.close();
923  gridView_.comm().barrier();
924 
925  // if we are rank 0, write .pvtu/.pvtp parallel header
926  fullname = getParallelHeaderName(name, path, commSize);
927  if( commRank ==0 )
928  {
929  file.open(fullname.c_str());
930  if (! file.is_open())
931  DUNE_THROW(IOError, "Could not write to parallel file " << fullname);
932  writeParallelHeader(file,name,relpiecepath, commSize );
933  file.close();
934  }
935  gridView_.comm().barrier();
936  return fullname;
937  }
938 
939  private:
941 
958  void writeParallelHeader(std::ostream& s, const std::string& piecename,
959  const std::string& piecepath, const int commSize)
960  {
961  VTK::FileType fileType =
963 
964  VTK::PVTUWriter writer(s, fileType);
965 
966  writer.beginMain();
967 
968  // PPointData
969  {
970  std::string scalars, vectors;
971  std::tie(scalars,vectors) = getDataNames(vertexdata);
972  writer.beginPointData(scalars, vectors);
973  }
974  for (auto it = vertexdata.begin(),
975  end = vertexdata.end();
976  it != end;
977  ++it)
978  {
979  unsigned writecomps = it->fieldInfo().size();
980  if(writecomps == 2) writecomps = 3;
981  writer.addArray<float>(it->name(), writecomps);
982  }
983  writer.endPointData();
984 
985  // PCellData
986  {
987  std::string scalars, vectors;
988  std::tie(scalars,vectors) = getDataNames(celldata);
989  writer.beginCellData(scalars, vectors);
990  }
991  for (auto it = celldata.begin(),
992  end = celldata.end();
993  it != end;
994  ++it)
995  {
996  unsigned writecomps = it->fieldInfo().size();
997  if(writecomps == 2) writecomps = 3;
998  writer.addArray<float>(it->name(), writecomps);
999  }
1000  writer.endCellData();
1001 
1002  // PPoints
1003  writer.beginPoints();
1004  writer.addArray<float>("Coordinates", 3);
1005  writer.endPoints();
1006 
1007  // Pieces
1008  for( int i = 0; i < commSize; ++i )
1009  {
1010  const std::string& fullname = getParallelPieceName(piecename,
1011  piecepath, i,
1012  commSize);
1013  writer.addPiece(fullname);
1014  }
1015 
1016  writer.endMain();
1017  }
1018 
1020  void writeDataFile (std::ostream& s)
1021  {
1022  VTK::FileType fileType =
1023  (n == 1) ? VTK::polyData : VTK::unstructuredGrid;
1024 
1025  VTK::VTUWriter writer(s, outputtype, fileType);
1026 
1027  // Grid characteristics
1028  vertexmapper = new VertexMapper( gridView_ );
1029  if (datamode == VTK::conforming)
1030  {
1031  number.resize(vertexmapper->size());
1032  for (std::vector<int>::size_type i=0; i<number.size(); i++) number[i] = -1;
1033  }
1035 
1036  writer.beginMain(ncells, nvertices);
1037  writeAllData(writer);
1038  writer.endMain();
1039 
1040  // write appended binary data section
1041  if(writer.beginAppended())
1042  writeAllData(writer);
1043  writer.endAppended();
1044 
1045  delete vertexmapper; number.clear();
1046  }
1047 
1048  void writeAllData(VTK::VTUWriter& writer) {
1049  // PointData
1050  writeVertexData(writer);
1051 
1052  // CellData
1053  writeCellData(writer);
1054 
1055  // Points
1056  writeGridPoints(writer);
1057 
1058  // Cells
1059  writeGridCells(writer);
1060  }
1061 
1062  protected:
1063  std::string getFormatString() const
1064  {
1065  if (outputtype==VTK::ascii)
1066  return "ascii";
1067  if (outputtype==VTK::base64)
1068  return "binary";
1070  return "appended";
1072  return "appended";
1073  DUNE_THROW(IOError, "VTKWriter: unsupported OutputType" << outputtype);
1074  }
1075 
1076  std::string getTypeString() const
1077  {
1078  if (n==1)
1079  return "PolyData";
1080  else
1081  return "UnstructuredGrid";
1082  }
1083 
1085  virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
1086  {
1087  nvertices = 0;
1088  ncells = 0;
1089  ncorners = 0;
1090  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1091  {
1092  ncells++;
1093  // because of the use of vertexmapper->map(), this iteration must be
1094  // in the order of Dune's numbering.
1095  const int subEntities = it->subEntities(n);
1096  for (int i=0; i<subEntities; ++i)
1097  {
1098  ncorners++;
1099  if (datamode == VTK::conforming)
1100  {
1101  int alpha = vertexmapper->subIndex(*it,i,n);
1102  if (number[alpha]<0)
1103  number[alpha] = nvertices++;
1104  }
1105  else
1106  {
1107  nvertices++;
1108  }
1109  }
1110  }
1111  }
1112 
1113  template<typename T>
1114  std::tuple<std::string,std::string> getDataNames(const T& data) const
1115  {
1116  std::string scalars = "";
1117  for (auto it = data.begin(),
1118  end = data.end();
1119  it != end;
1120  ++it)
1121  if (it->fieldInfo().type() == VTK::FieldInfo::Type::scalar)
1122  {
1123  scalars = it->name();
1124  break;
1125  }
1126 
1127  std::string vectors = "";
1128  for (auto it = data.begin(),
1129  end = data.end();
1130  it != end;
1131  ++it)
1132  if (it->fieldInfo().type() == VTK::FieldInfo::Type::vector)
1133  {
1134  vectors = it->name();
1135  break;
1136  }
1137  return std::make_tuple(scalars,vectors);
1138  }
1139 
1140  template<typename Data, typename Iterator>
1141  void writeData(VTK::VTUWriter& writer, const Data& data, const Iterator begin, const Iterator end, int nentries)
1142  {
1143  for (auto it = data.begin(),
1144  iend = data.end();
1145  it != iend;
1146  ++it)
1147  {
1148  const auto& f = *it;
1149  VTK::FieldInfo fieldInfo = f.fieldInfo();
1150  std::size_t writecomps = fieldInfo.size();
1151  switch (fieldInfo.type())
1152  {
1154  break;
1156  // vtk file format: a vector data always should have 3 comps (with
1157  // 3rd comp = 0 in 2D case)
1158  if (writecomps > 3)
1159  DUNE_THROW(IOError,"Cannot write VTK vectors with more than 3 components (components was " << writecomps << ")");
1160  writecomps = 3;
1161  break;
1163  DUNE_THROW(NotImplemented,"VTK output for tensors not implemented yet");
1164  }
1165  std::shared_ptr<VTK::DataArrayWriter<float> > p
1166  (writer.makeArrayWriter<float>(f.name(), writecomps, nentries));
1167  if(!p->writeIsNoop())
1168  for (Iterator eit = begin; eit!=end; ++eit)
1169  {
1170  const Entity & e = *eit;
1171  f.bind(e);
1172  f.write(eit.position(),*p);
1173  f.unbind();
1174  // vtk file format: a vector data always should have 3 comps
1175  // (with 3rd comp = 0 in 2D case)
1176  for (std::size_t j=fieldInfo.size(); j < writecomps; ++j)
1177  p->write(0.0);
1178  }
1179  }
1180  }
1181 
1183  virtual void writeCellData(VTK::VTUWriter& writer)
1184  {
1185  if(celldata.size() == 0)
1186  return;
1187 
1188  std::string scalars, vectors;
1189  std::tie(scalars,vectors) = getDataNames(celldata);
1190 
1191  writer.beginCellData(scalars, vectors);
1193  writer.endCellData();
1194  }
1195 
1197  virtual void writeVertexData(VTK::VTUWriter& writer)
1198  {
1199  if(vertexdata.size() == 0)
1200  return;
1201 
1202  std::string scalars, vectors;
1203  std::tie(scalars,vectors) = getDataNames(vertexdata);
1204 
1205  writer.beginPointData(scalars, vectors);
1207  writer.endPointData();
1208  }
1209 
1211  virtual void writeGridPoints(VTK::VTUWriter& writer)
1212  {
1213  writer.beginPoints();
1214 
1215  std::shared_ptr<VTK::DataArrayWriter<float> > p
1216  (writer.makeArrayWriter<float>("Coordinates", 3, nvertices));
1217  if(!p->writeIsNoop()) {
1218  VertexIterator vEnd = vertexEnd();
1219  for (VertexIterator vit=vertexBegin(); vit!=vEnd; ++vit)
1220  {
1221  int dimw=w;
1222  for (int j=0; j<std::min(dimw,3); j++)
1223  p->write((*vit).geometry().corner(vit.localindex())[j]);
1224  for (int j=std::min(dimw,3); j<3; j++)
1225  p->write(0.0);
1226  }
1227  }
1228  // free the VTK::DataArrayWriter before touching the stream
1229  p.reset();
1230 
1231  writer.endPoints();
1232  }
1233 
1235  virtual void writeGridCells(VTK::VTUWriter& writer)
1236  {
1237  writer.beginCells();
1238 
1239  // connectivity
1240  {
1241  std::shared_ptr<VTK::DataArrayWriter<int> > p1
1242  (writer.makeArrayWriter<int>("connectivity", 1, ncorners));
1243  if(!p1->writeIsNoop())
1244  for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
1245  p1->write(it.id());
1246  }
1247 
1248  // offsets
1249  {
1250  std::shared_ptr<VTK::DataArrayWriter<int> > p2
1251  (writer.makeArrayWriter<int>("offsets", 1, ncells));
1252  if(!p2->writeIsNoop()) {
1253  int offset = 0;
1254  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1255  {
1256  offset += it->subEntities(n);
1257  p2->write(offset);
1258  }
1259  }
1260  }
1261 
1262  // types
1263  if (n>1)
1264  {
1265  std::shared_ptr<VTK::DataArrayWriter<unsigned char> > p3
1266  (writer.makeArrayWriter<unsigned char>("types", 1, ncells));
1267  if(!p3->writeIsNoop())
1268  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1269  {
1270  int vtktype = VTK::geometryType(it->type());
1271  p3->write(vtktype);
1272  }
1273  }
1274 
1275  writer.endCells();
1276  }
1277 
1278  protected:
1279  // the list of registered functions
1280  std::list<VTKLocalFunction> celldata;
1281  std::list<VTKLocalFunction> vertexdata;
1282 
1283  // the grid
1285 
1286  // temporary grid information
1287  int ncells;
1290  private:
1291  VertexMapper* vertexmapper;
1292  // in conforming mode, for each vertex id (as obtained by vertexmapper)
1293  // hold its number in the iteration order (VertexIterator)
1294  std::vector<int> number;
1295  VTK::DataMode datamode;
1296  protected:
1298  };
1299 
1300 }
1301 
1302 #endif
Dune::P0VTKFunction
Take a vector and interpret it as cell data for the VTKWriter.
Definition: function.hh:88
Dune::VTK::DataArrayWriter::write
virtual void write(T data)=0
write one data element
Dune::Alberta::min
int min(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:346
Dune::VTKWriter::VTKLocalFunction::VTKFunctionWrapper::bind
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:231
Dune::VTKWriter::VTKLocalFunction::VTKFunctionWrapper::VTKFunctionWrapper
VTKFunctionWrapper(const std::shared_ptr< const VTKFunction > &f)
Definition: vtkwriter.hh:226
pvtuwriter.hh
Dune::VTK::VTUWriter::endCells
void endCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:283
Dune::VTKWriter::vertexBegin
VertexIterator vertexBegin() const
Definition: vtkwriter.hh:449
Dune::VTK::VTUWriter::endCellData
void endCellData()
finish CellData section
Definition: vtuwriter.hh:218
Dune::VTKWriter::celldata
std::list< VTKLocalFunction > celldata
Definition: vtkwriter.hh:1280
Dune::VTKWriter::getDataNames
std::tuple< std::string, std::string > getDataNames(const T &data) const
Definition: vtkwriter.hh:1114
vtuwriter.hh
Dune::VTKWriter::VertexIterator::basicIncrement
void basicIncrement()
Definition: vtkwriter.hh:381
Dune
Include standard header files.
Definition: agrid.hh:59
Dune::VTKWriter::pwrite
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:893
Dune::VTK::base64
Output to the file is inline base64 binary.
Definition: common.hh:44
Dune::VTKWriter::addVertexData
void addVertexData(F &&f, VTK::FieldInfo vtkFieldInfo)
Definition: vtkwriter.hh:636
Dune::VTKWriter::CornerIterator
Iterate over the elements' corners.
Definition: vtkwriter.hh:478
Dune::VTKWriter::VTKLocalFunction::FunctionWrapperBase::bind
virtual void bind(const Entity &e)=0
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Dune::VTKWriter::VTKLocalFunction::VTKLocalFunction
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:256
Dune::VTK::VTUWriter::beginCellData
void beginCellData(const std::string &scalars="", const std::string &vectors="")
start CellData section
Definition: vtuwriter.hh:203
Dune::VTK::FieldInfo::Type::tensor
tensor field (always 3x3)
Dune::VTK::OutputType
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:40
Dune::VTK::VTUWriter::endPoints
void endPoints()
finish section for the point coordinates
Definition: vtuwriter.hh:247
Dune::VTKWriter::write
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:833
Dune::VTKWriter::countEntities
virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
count the vertices, cells and corners
Definition: vtkwriter.hh:1085
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper::unbind
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:191
Dune::VTKWriter::VertexIterator::equals
bool equals(const VertexIterator &cit) const
Definition: vtkwriter.hh:426
Dune::VTK::appendedraw
Output is to the file is appended raw binary.
Definition: common.hh:46
Dune::VTKWriter::addVertexData
void addVertexData(const Container &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:658
Dune::VTK::polyData
for .vtp files (PolyData)
Definition: common.hh:292
Dune::VTKWriter::CornerIterator::equals
bool equals(const CornerIterator &cit) const
Definition: vtkwriter.hh:523
Dune::VTK::conforming
Output conforming data.
Definition: common.hh:70
Dune::VTKWriter::VTKLocalFunction::name
std::string name() const
Returns the name of the data set.
Definition: vtkwriter.hh:283
Dune::VTKWriter::VTKLocalFunction::unbind
void unbind() const
Unbind the data set from the currently bound entity.
Definition: vtkwriter.hh:301
Dune::VTKSequenceWriter
Writer for the ouput of grid functions in the vtk format.
Definition: vtksequencewriter.hh:26
Dune::VTKWriter::vertexEnd
VertexIterator vertexEnd() const
Definition: vtkwriter.hh:456
Dune::VTK::FieldInfo::Type::vector
vector-valued field (always 3D, will be padded if necessary)
Dune::VTKWriter::writeGridPoints
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: vtkwriter.hh:1211
Dune::GridView::comm
const CollectiveCommunication & comm() const
obtain collective communication object
Definition: common/gridview.hh:247
Dune::VTKWriter::VertexIterator::position
const FieldVector< DT, n > & position() const
position of vertex inside the entity
Definition: vtkwriter.hh:442
Dune::VTKWriter::VTKWriter
VTKWriter(const GridView &gridView, VTK::DataMode dm=VTK::conforming)
Construct a VTKWriter working on a specific GridView.
Definition: vtkwriter.hh:576
Dune::MultipleCodimMultipleGeomTypeMapper::subIndex
Index subIndex(const typename GV::template Codim< 0 >::Entity &e, int i, unsigned int codim) const
Map subentity of codim 0 entity to array index.
Definition: mcmgmapper.hh:160
mcmgmapper.hh
Mapper for multiple codim and multiple geometry types.
Dune::VTKWriter::addVertexData
void addVertexData(const std::shared_ptr< const VTKFunction > &p)
Add a grid function that lives on the vertices of the grid to the visualization.
Definition: vtkwriter.hh:630
Dune::VTKWriter::addCellData
void addCellData(const Container &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:613
Dune::VTK::VTUWriter::endPointData
void endPointData()
finish PointData section
Definition: vtuwriter.hh:180
Dune::VTKWriter::VTKLocalFunction::fieldInfo
const VTK::FieldInfo & fieldInfo() const
Returns the VTK::FieldInfo for the data set.
Definition: vtkwriter.hh:289
Dune::PartitionType
PartitionType
Attributes used in the generic overlap model.
Definition: gridenums.hh:28
Dune::VTKWriter::gridView_
GridView gridView_
Definition: vtkwriter.hh:1284
Dune::VTK::DataArrayWriter< float >
Dune::VTK::FieldInfo::name
std::string name() const
The name of the data field.
Definition: common.hh:330
gridenums.hh
Dune::VTK::FieldInfo
Descriptor struct for VTK fields.
Definition: common.hh:306
streams.hh
Dune::VTKWriter::cornerBegin
CornerIterator cornerBegin() const
Definition: vtkwriter.hh:554
Dune::VTKWriter::CornerIterator::dereference
EntityReference dereference() const
Definition: vtkwriter.hh:529
Dune::VTK::FieldInfo::Type::scalar
Dune::VTK::ascii
Output to the file is in ascii.
Definition: common.hh:42
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper::bind
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:186
Dune::GridView::IndexSet
Traits ::IndexSet IndexSet
type of the index set
Definition: common/gridview.hh:81
Dune::VTK::VTUWriter::beginPoints
void beginPoints()
start section for the point coordinates
Definition: vtuwriter.hh:236
Dune::PartitionIteratorType
PartitionIteratorType
Parameter to be used for the parallel level- and leaf iterators.
Definition: gridenums.hh:134
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper::write
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:196
Dune::VTKWriter::CornerIterator::CornerIterator
CornerIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm, const std::vector< int > &num)
Definition: vtkwriter.hh:499
Dune::VTK::nonconforming
Output non-conforming data.
Definition: common.hh:78
Dune::VTK::geometryType
GeometryType geometryType(const Dune::GeometryType &t)
mapping from GeometryType to VTKGeometryType
Definition: common.hh:195
Dune::VTKWriter::clear
void clear()
clear list of registered functions
Definition: vtkwriter.hh:672
Dune::VTKWriter::write
std::string write(const std::string &name, VTK::OutputType type=VTK::ascii)
write output (interface might change later)
Definition: vtkwriter.hh:695
Dune::VTKWriter::VTKLocalFunction::FunctionWrapperBase::~FunctionWrapperBase
virtual ~FunctionWrapperBase()
Definition: vtkwriter.hh:169
Dune::GridView::dimension
The dimension of the grid.
Definition: common/gridview.hh:128
Dune::MultipleCodimMultipleGeomTypeMapper
Implementation class for a multiple codim and multiple geometry type mapper.
Definition: mcmgmapper.hh:101
Dune::VTK::VTUWriter::beginPointData
void beginPointData(const std::string &scalars="", const std::string &vectors="")
start PointData section
Definition: vtuwriter.hh:165
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper::Function
typename std::decay< F >::type Function
Definition: vtkwriter.hh:179
Dune::VTKWriter::writeCellData
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: vtkwriter.hh:1183
Dune::VTK::VTUWriter::beginCells
void beginCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:272
Dune::VTKWriter::VertexIterator::VertexIterator
VertexIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm)
Definition: vtkwriter.hh:398
Dune::VTKWriter::VTKLocalFunction::_fieldInfo
VTK::FieldInfo _fieldInfo
Definition: vtkwriter.hh:313
Dune::VTKWriter::CornerIterator::id
int id() const
Process-local consecutive zero-starting vertex id.
Definition: vtkwriter.hh:538
Dune::VTK::appendedbase64
Output is to the file is appended base64 binary.
Definition: common.hh:48
Dune::VTKWriter::VTKLocalFunction::FunctionWrapperBase
Base class for polymorphic container of underlying data set.
Definition: vtkwriter.hh:154
Dune::VTKWriter::cellBegin
CellIterator cellBegin() const
Definition: vtkwriter.hh:338
Dune::VTKWriter::addCellData
void addCellData(const std::shared_ptr< const VTKFunction > &p)
Add a grid function that lives on the cells of the grid to the visualization.
Definition: vtkwriter.hh:586
Dune::VTKWriter::vertexdata
std::list< VTKLocalFunction > vertexdata
Definition: vtkwriter.hh:1281
Dune::VTKWriter::getParallelHeaderName
std::string getParallelHeaderName(const std::string &name, const std::string &path, int commSize) const
return name of a parallel header file
Definition: vtkwriter.hh:777
Dune::VTK::VTUWriter
Dump a .vtu/.vtp files contents to a stream.
Definition: vtuwriter.hh:96
Dune::VTK::PVTUWriter
Dump a .vtu/.vtp files contents to a stream.
Definition: pvtuwriter.hh:60
Dune::All_Partition
all entities
Definition: gridenums.hh:139
Dune::VTKWriter::VTKLocalFunction::bind
void bind(const Entity &e) const
Bind the data set to grid entity e.
Definition: vtkwriter.hh:295
Dune::VTKWriter::CellIterator::CellIterator
CellIterator(const GridCellIterator &x)
construct a CellIterator from the gridview's Iterator.
Definition: vtkwriter.hh:329
Dune::GridView::dimensionworld
The dimension of the world the grid lives in.
Definition: common/gridview.hh:132
Dune::VTK::renumber
int renumber(const Dune::GeometryType &t, int i)
renumber VTK <-> Dune
Definition: common.hh:224
Dune::VTK::DataMode
DataMode
Whether to produce conforming or non-conforming output.
Definition: common.hh:64
Dune::VTKWriter::VTKLocalFunction::VTKLocalFunction
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:264
Dune::GridView::Grid
Traits ::Grid Grid
type of the grid
Definition: common/gridview.hh:78
Dune::VTKFunction
A base class for grid functions with any return type and dimension.
Definition: function.hh:38
Dune::detail::_has_local_context
Definition: vtkwriter.hh:53
Dune::VTKWriter::VTKLocalFunction::VTKLocalFunction
VTKLocalFunction(const std::shared_ptr< const VTKFunction > &vtkFunctionPtr)
Construct a VTKLocalFunction for a legacy VTKFunction.
Definition: vtkwriter.hh:273
Dune::VTKWriter::~VTKWriter
virtual ~VTKWriter()
destructor
Definition: vtkwriter.hh:679
Dune::VTKWriter::VTKLocalFunction::_f
std::shared_ptr< FunctionWrapperBase > _f
Definition: vtkwriter.hh:312
common.hh
Common stuff for the VTKWriter.
Dune::VTKWriter::VTKLocalFunction::FunctionWrapperBase::write
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.
Dune::InteriorBorder_Partition
interior and border entities
Definition: gridenums.hh:136
Dune::VTK::VTUWriter::makeArrayWriter
DataArrayWriter< T > * makeArrayWriter(const std::string &name, unsigned ncomps, unsigned nitems)
acquire a DataArrayWriter
Definition: vtuwriter.hh:379
Dune::InteriorEntity
all interior entities
Definition: gridenums.hh:29
Dune::VTKWriter::VertexIterator
Iterate over the grid's vertices.
Definition: vtkwriter.hh:363
Dune::VTKWriter::CellIterator::position
const FieldVector< DT, n > position() const
Definition: vtkwriter.hh:332
Dune::VTKWriter::VTKLocalFunction::FunctionWrapperBase::unbind
virtual void unbind()=0
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Dune::VTKWriter::getParallelPieceName
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:746
Dune::VTKWriter::CornerIterator::increment
void increment()
Definition: vtkwriter.hh:507
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper::FunctionWrapper
FunctionWrapper(F_ &&f)
Definition: vtkwriter.hh:182
Dune::VTKWriter::getFormatString
std::string getFormatString() const
Definition: vtkwriter.hh:1063
Dune::VTKWriter::VertexIterator::increment
void increment()
Definition: vtkwriter.hh:409
Dune::VTKWriter::getSerialPieceName
std::string getSerialPieceName(const std::string &name, const std::string &path) const
return name of a serial piece file
Definition: vtkwriter.hh:809
Dune::VTKWriter::VTKLocalFunction::VTKFunctionWrapper
Type erasure implementation for legacy VTKFunctions.
Definition: vtkwriter.hh:223
dataarraywriter.hh
Data array writers for the VTKWriter.
Dune::VTKWriter::VTKFunction
Dune::VTKFunction< GridView > VTKFunction
Definition: vtkwriter.hh:138
Dune::VTKWriter::FunctionIterator
std::list< VTKLocalFunction >::const_iterator FunctionIterator
Definition: vtkwriter.hh:317
Dune::VTKWriter::VTKLocalFunction
Type erasure wrapper for VTK data sets.
Definition: vtkwriter.hh:146
Dune::MultipleCodimMultipleGeomTypeMapper::size
int size() const
Return total number of entities in the entity set managed by the mapper.
Definition: mcmgmapper.hh:175
Dune::VTKWriter::nvertices
int nvertices
Definition: vtkwriter.hh:1288
Dune::VTKWriter::VertexIterator::dereference
EntityReference dereference() const
Definition: vtkwriter.hh:432
Dune::VTK::FieldInfo::type
Type type() const
The type of the data field.
Definition: common.hh:336
Dune::GridView::ctype
Grid::ctype ctype
type used for coordinates in grid
Definition: common/gridview.hh:125
Dune::VTKWriter::CellIterator
Iterator over the grids elements.
Definition: vtkwriter.hh:325
Dune::VTKWriter::VTKLocalFunction::write
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:307
Dune::VTKWriter::VertexIterator::localindex
int localindex() const
index of vertex within the entity, in Dune-numbering
Definition: vtkwriter.hh:437
Dune::VTK::unstructuredGrid
for .vtu files (UnstructuredGrid)
Definition: common.hh:294
Dune::VTKWriter::getTypeString
std::string getTypeString() const
Definition: vtkwriter.hh:1076
Dune::P1VTKFunction
Take a vector and interpret it as point data for the VTKWriter.
Definition: function.hh:186
Dune::VTKWriter::writeVertexData
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: vtkwriter.hh:1197
Dune::VTKWriteTypeTraits::IsLocalFunction
Definition: vtkwriter.hh:66
Dune::VTKWriter::VTKLocalFunction::Writer
VTK::DataArrayWriter< float > Writer
Definition: vtkwriter.hh:151
Dune::VTK::FieldInfo::size
std::size_t size() const
The number of components in the data field.
Definition: common.hh:342
Dune::VTKWriter::cellEnd
CellIterator cellEnd() const
Definition: vtkwriter.hh:343
Dune::VTKWriter::addCellData
void addCellData(F &&f, VTK::FieldInfo vtkFieldInfo)
Definition: vtkwriter.hh:592
Dune::VTKSequenceWriterBase
Base class to write pvd-files which contains a list of all collected vtk-files.
Definition: vtksequencewriterbase.hh:31
Dune::VTK::FileType
FileType
which type of VTK file to write
Definition: common.hh:290
Dune::VTKWriter::VTKLocalFunction::VTKFunctionWrapper::unbind
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:236
Dune::VTKWriter::VTKLocalFunction::FunctionWrapper
Type erasure implementation for functions conforming to the dune-functions LocalFunction interface.
Definition: vtkwriter.hh:176
Dune::VTKWriter::pwrite
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:727
Dune::VTKWriter
Writer for the ouput of grid functions in the vtk format.
Definition: vtkwriter.hh:86
Dune::GridView
Grid view abstract base class.
Definition: common/gridview.hh:59
Dune::VTKWriter::writeData
void writeData(VTK::VTUWriter &writer, const Data &data, const Iterator begin, const Iterator end, int nentries)
Definition: vtkwriter.hh:1141
Dune::VTKWriter::ncorners
int ncorners
Definition: vtkwriter.hh:1289
Dune::VTKWriter::cornerEnd
CornerIterator cornerEnd() const
Definition: vtkwriter.hh:561
Dune::VTKWriter::VTKLocalFunction::VTKFunctionWrapper::write
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:241
function.hh
Functions for VTK output.
Dune::VTKWriter::ncells
int ncells
Definition: vtkwriter.hh:1287
Dune::VTKWriter::outputtype
VTK::OutputType outputtype
Definition: vtkwriter.hh:1297
Dune::VTKWriter::writeGridCells
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: vtkwriter.hh:1235