dune-grid  2.7.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 #include <memory>
14 
15 #include <vector>
16 #include <list>
17 #include <map>
18 
19 #include <dune/common/typetraits.hh>
20 #include <dune/common/exceptions.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 Impl
51  {
52  // Check whether type F has a method 'bind' (see the dune-functions interface)
53  template< class F, class E, class = void >
54  struct IsBindable
55  : std::false_type
56  {};
57 
58  template< class F, class E >
59  struct IsBindable< F, E, void_t< decltype( std::declval< F & >().bind( std::declval< const E & >() ) ), decltype( std::declval< F & >().unbind() ) > >
60  : std::true_type
61  {};
62 
63  // Check whether localFunction(F) can be called (see the dune-functions interface)
64  template< class F, class = void >
65  struct HasLocalFunction
66  : std::false_type
67  {};
68 
69  template< class F >
70  struct HasLocalFunction< F, void_t< decltype( localFunction( std::declval< F& >() ) ) > >
71  : std::true_type
72  {};
73 
74  } // namespace Impl
75 
76  // Forward-declaration here, so the class can be friend of VTKWriter
77  template <class GridView>
79  template <class GridView>
80  class VTKSequenceWriter;
81 
90  template< class GridView >
91  class VTKWriter {
92 
93  // VTKSequenceWriterBase needs getSerialPieceName
94  // and getParallelHeaderName
95  friend class VTKSequenceWriterBase<GridView>;
96  // VTKSequenceWriter needs the grid view, to get the MPI size and rank
97  friend class VTKSequenceWriter<GridView>;
98 
99  // extract types
100  typedef typename GridView::Grid Grid;
101  typedef typename GridView::ctype DT;
102  enum { n = GridView::dimension };
103  enum { w = GridView::dimensionworld };
104 
105  typedef typename GridView::template Codim< 0 >::Entity Cell;
106  typedef typename GridView::template Codim< n >::Entity Vertex;
107  typedef Cell Entity;
108 
109  typedef typename GridView::IndexSet IndexSet;
110 
111  static const PartitionIteratorType VTK_Partition = InteriorBorder_Partition;
112  //static const PartitionIteratorType VTK_Partition = All_Partition;
113 
114  typedef typename GridView::template Codim< 0 >
115  ::template Partition< VTK_Partition >::Iterator
116  GridCellIterator;
117  typedef typename GridView::template Codim< n >
118  ::template Partition< VTK_Partition >::Iterator
119  GridVertexIterator;
120 
121  typedef typename GridCellIterator::Reference EntityReference;
122 
123  typedef typename GridView::template Codim< 0 >
124  ::Entity::Geometry::LocalCoordinate Coordinate;
125 
127 
128  // return true if entity should be skipped in Vertex and Corner iterator
129  static bool skipEntity( const PartitionType entityType )
130  {
131  switch( VTK_Partition )
132  {
133  // for All_Partition no entity has to be skipped
134  case All_Partition: return false;
135  case InteriorBorder_Partition: return ( entityType != InteriorEntity );
136  default: DUNE_THROW(NotImplemented,"Add check for this partition type");
137  }
138  return false ;
139  }
140 
141  public:
142 
144 
145  protected:
146 
148 
152  {
153 
154  public:
155 
157 
160  {
161 
163  virtual void bind(const Entity& e) = 0;
164 
166  virtual void unbind() = 0;
167 
169 
172  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const = 0;
173 
175  {}
176 
177  };
178 
180  template<typename F>
182  : public FunctionWrapperBase
183  {
184  using Function = typename std::decay<F>::type;
185 
186  template<typename F_>
188  : _f(std::forward<F_>(f))
189  {}
190 
191  virtual void bind(const Entity& e)
192  {
193  _f.bind(e);
194  }
195 
196  virtual void unbind()
197  {
198  _f.unbind();
199  }
200 
201  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
202  {
203  auto r = _f(pos);
204  // we need to do different things here depending on whether r supports indexing into it or not.
205  do_write(w,r,count,IsIndexable<decltype(r)>());
206  }
207 
208  private:
209 
210  template<typename R>
211  void do_write(Writer& w, const R& r, std::size_t count, std::true_type) const
212  {
213  for (std::size_t i = 0; i < count; ++i)
214  w.write(r[i]);
215  }
216 
217  template<typename R>
218  void do_write(Writer& w, const R& r, std::size_t count, std::false_type) const
219  {
220  assert(count == 1);
221  w.write(r);
222  }
223 
224  Function _f;
225  };
226 
228  template<typename F>
230  : public FunctionWrapperBase
231  {
232  using Function = typename std::decay<F>::type;
233 
234  template<typename F_>
236  : _f(std::forward<F_>(f))
237  , element_(nullptr)
238  {}
239 
240  virtual void bind(const Entity& e)
241  {
242  element_ = &e;
243  }
244 
245  virtual void unbind()
246  {
247  element_ = nullptr;
248  }
249 
250  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
251  {
252  auto globalPos = element_->geometry().global(pos);
253  auto r = _f(globalPos);
254  Hybrid::ifElse(IsIndexable<decltype(r)>(),
255  [&](auto id) {
256  for (std::size_t i = 0; i < count; ++i)
257  w.write(id(r)[i]);
258  },
259  [&](auto id) {
260  assert(count == 1);
261  w.write(id(r));
262  });
263  }
264  private:
265  Function _f;
266  const Entity* element_;
267  };
268 
271  : public FunctionWrapperBase
272  {
273  VTKFunctionWrapper(const std::shared_ptr< const VTKFunction >& f)
274  : _f(f)
275  , _entity(nullptr)
276  {}
277 
278  virtual void bind(const Entity& e)
279  {
280  _entity = &e;
281  }
282 
283  virtual void unbind()
284  {
285  _entity = nullptr;
286  }
287 
288  virtual void write(const Coordinate& pos, Writer& w, std::size_t count) const
289  {
290  for (std::size_t i = 0; i < count; ++i)
291  w.write(_f->evaluate(i,*_entity,pos));
292  }
293 
294  private:
295 
296  std::shared_ptr< const VTKFunction > _f;
297  const Entity* _entity;
298 
299  };
300 
302  template<typename F, std::enable_if_t<Impl::IsBindable<F, Entity>::value, int> = 0>
304  : _f(std::make_unique<FunctionWrapper<F> >(std::forward<F>(f)))
306  {}
307 
309  // That is, a function that you can create a LocalFunction for, and evaluate that in element coordinates
310  template<typename F, std::enable_if_t<not Impl::IsBindable<F, Entity>::value && Impl::HasLocalFunction<F>::value, int> = 0>
312  : _f(std::make_unique< FunctionWrapper<
313  typename std::decay<decltype(localFunction(std::forward<F>(f)))>::type
314  > >(localFunction(std::forward<F>(f))))
316  {}
317 
319  // That is, a function that can be evaluated in global coordinates of the domain
320  template<typename F, std::enable_if_t<not Impl::IsBindable<F, Entity>::value && not Impl::HasLocalFunction<F>::value, int> = 0>
322  : _f(std::make_unique< GlobalFunctionWrapper<F> >(std::forward<F>(f)))
324  {}
325 
327  explicit VTKLocalFunction (const std::shared_ptr< const VTKFunction >& vtkFunctionPtr)
328  : _f(std::make_unique<VTKFunctionWrapper>(vtkFunctionPtr))
329  , _fieldInfo(
330  vtkFunctionPtr->name(),
331  vtkFunctionPtr->ncomps() > 1 ? VTK::FieldInfo::Type::vector : VTK::FieldInfo::Type::scalar,
332  vtkFunctionPtr->ncomps(),
333  vtkFunctionPtr->precision()
334  )
335  {}
336 
338  std::string name() const
339  {
340  return fieldInfo().name();
341  }
342 
344  const VTK::FieldInfo& fieldInfo() const
345  {
346  return _fieldInfo;
347  }
348 
350  void bind(const Entity& e) const
351  {
352  _f->bind(e);
353  }
354 
356  void unbind() const
357  {
358  _f->unbind();
359  }
360 
362  void write(const Coordinate& pos, Writer& w) const
363  {
364  _f->write(pos,w,fieldInfo().size());
365  }
366 
367  std::shared_ptr<FunctionWrapperBase> _f;
369 
370  };
371 
372  typedef typename std::list<VTKLocalFunction>::const_iterator FunctionIterator;
373 
375 
380  class CellIterator : public GridCellIterator
381  {
382  public:
384  CellIterator(const GridCellIterator & x) : GridCellIterator(x) {}
387  const FieldVector<DT,n> position() const
388  {
389  return ReferenceElements<DT,n>::general((*this)->type()).position(0,0);
390  }
391  };
392 
394  {
395  return gridView_.template begin< 0, VTK_Partition >();
396  }
397 
399  {
400  return gridView_.template end< 0, VTK_Partition >();
401  }
402 
404 
419  public ForwardIteratorFacade<VertexIterator, const Entity, EntityReference, int>
420  {
421  GridCellIterator git;
422  GridCellIterator gend;
423  VTK::DataMode datamode;
424  // Index of the currently visited corner within the current element.
425  // NOTE: this is in Dune-numbering, in contrast to CornerIterator.
426  int cornerIndexDune;
427  const VertexMapper & vertexmapper;
428  std::vector<bool> visited;
429  // in conforming mode, for each vertex id (as obtained by vertexmapper)
430  // hold its number in the iteration order (VertexIterator)
431  int offset;
432 
433  // hide operator ->
434  void operator->();
435  protected:
437  {
438  if( git == gend )
439  return;
440  ++cornerIndexDune;
441  const int numCorners = git->subEntities(n);
442  if( cornerIndexDune == numCorners )
443  {
444  offset += numCorners;
445  cornerIndexDune = 0;
446 
447  ++git;
448  while( (git != gend) && skipEntity( git->partitionType() ) )
449  ++git;
450  }
451  }
452  public:
453  VertexIterator(const GridCellIterator & x,
454  const GridCellIterator & end,
455  const VTK::DataMode & dm,
456  const VertexMapper & vm) :
457  git(x), gend(end), datamode(dm), cornerIndexDune(0),
458  vertexmapper(vm), visited(vm.size(), false),
459  offset(0)
460  {
461  if (datamode == VTK::conforming && git != gend)
462  visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
463  }
464  void increment ()
465  {
466  switch (datamode)
467  {
468  case VTK::conforming :
469  while(visited[vertexmapper.subIndex(*git,cornerIndexDune,n)])
470  {
471  basicIncrement();
472  if (git == gend) return;
473  }
474  visited[vertexmapper.subIndex(*git,cornerIndexDune,n)] = true;
475  break;
476  case VTK::nonconforming :
477  basicIncrement();
478  break;
479  }
480  }
481  bool equals (const VertexIterator & cit) const
482  {
483  return git == cit.git
484  && cornerIndexDune == cit.cornerIndexDune
485  && datamode == cit.datamode;
486  }
487  EntityReference dereference() const
488  {
489  return *git;
490  }
492  int localindex () const
493  {
494  return cornerIndexDune;
495  }
497  FieldVector<DT,n> position () const
498  {
499  return referenceElement<DT,n>(git->type())
500  .position(cornerIndexDune,n);
501  }
502  };
503 
505  {
506  return VertexIterator( gridView_.template begin< 0, VTK_Partition >(),
507  gridView_.template end< 0, VTK_Partition >(),
508  datamode, *vertexmapper );
509  }
510 
512  {
513  return VertexIterator( gridView_.template end< 0, VTK_Partition >(),
514  gridView_.template end< 0, VTK_Partition >(),
515  datamode, *vertexmapper );
516  }
517 
519 
534  public ForwardIteratorFacade<CornerIterator, const Entity, EntityReference, int>
535  {
536  GridCellIterator git;
537  GridCellIterator gend;
538  VTK::DataMode datamode;
539  // Index of the currently visited corner within the current element.
540  // NOTE: this is in VTK-numbering, in contrast to VertexIterator.
541  int cornerIndexVTK;
542  const VertexMapper & vertexmapper;
543  // in conforming mode, for each vertex id (as obtained by vertexmapper)
544  // hold its number in the iteration order of VertexIterator (*not*
545  // CornerIterator)
546  const std::vector<int> & number;
547  // holds the number of corners of all the elements we have seen so far,
548  // excluding the current element
549  int offset;
550 
551  // hide operator ->
552  void operator->();
553  public:
554  CornerIterator(const GridCellIterator & x,
555  const GridCellIterator & end,
556  const VTK::DataMode & dm,
557  const VertexMapper & vm,
558  const std::vector<int> & num) :
559  git(x), gend(end), datamode(dm), cornerIndexVTK(0),
560  vertexmapper(vm),
561  number(num), offset(0) {}
562  void increment ()
563  {
564  if( git == gend )
565  return;
566  ++cornerIndexVTK;
567  const int numCorners = git->subEntities(n);
568  if( cornerIndexVTK == numCorners )
569  {
570  offset += numCorners;
571  cornerIndexVTK = 0;
572 
573  ++git;
574  while( (git != gend) && skipEntity( git->partitionType() ) )
575  ++git;
576  }
577  }
578  bool equals (const CornerIterator & cit) const
579  {
580  return git == cit.git
581  && cornerIndexVTK == cit.cornerIndexVTK
582  && datamode == cit.datamode;
583  }
584  EntityReference dereference() const
585  {
586  return *git;
587  }
589 
593  int id () const
594  {
595  switch (datamode)
596  {
597  case VTK::conforming :
598  return
599  number[vertexmapper.subIndex(*git,VTK::renumber(*git,cornerIndexVTK),
600  n)];
601  case VTK::nonconforming :
602  return offset + VTK::renumber(*git,cornerIndexVTK);
603  default :
604  DUNE_THROW(IOError,"VTKWriter: unsupported DataMode" << datamode);
605  }
606  }
607  };
608 
610  {
611  return CornerIterator( gridView_.template begin< 0, VTK_Partition >(),
612  gridView_.template end< 0, VTK_Partition >(),
613  datamode, *vertexmapper, number );
614  }
615 
617  {
618  return CornerIterator( gridView_.template end< 0, VTK_Partition >(),
619  gridView_.template end< 0, VTK_Partition >(),
620  datamode, *vertexmapper, number );
621  }
622 
623  public:
632  explicit VTKWriter ( const GridView &gridView,
635  : gridView_( gridView ),
636  datamode( dm ),
637  coordPrec (coordPrecision),
638  polyhedralCellsPresent_( checkForPolyhedralCells() )
639  { }
640 
645  void addCellData (const std::shared_ptr< const VTKFunction > & p)
646  {
647  celldata.push_back(VTKLocalFunction(p));
648  }
649 
669  template<typename F>
670  void addCellData(F&& f, VTK::FieldInfo vtkFieldInfo)
671  {
672  celldata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
673  }
674 
690  template<class Container>
691  void addCellData (const Container& v, const std::string &name, int ncomps = 1,
693  {
694  typedef P0VTKFunction<GridView, Container> Function;
695  for (int c=0; c<ncomps; ++c) {
696  std::stringstream compName;
697  compName << name;
698  if (ncomps>1)
699  compName << "[" << c << "]";
700  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c, prec);
701  addCellData(std::shared_ptr< const VTKFunction >(p));
702  }
703  }
704 
709  void addVertexData (const std::shared_ptr< const VTKFunction > & p)
710  {
711  vertexdata.push_back(VTKLocalFunction(p));
712  }
713 
733  template<typename F>
734  void addVertexData(F&& f, VTK::FieldInfo vtkFieldInfo)
735  {
736  vertexdata.push_back(VTKLocalFunction(std::forward<F>(f),vtkFieldInfo));
737  }
738 
739 
755  template<class Container>
756  void addVertexData (const Container& v, const std::string &name, int ncomps=1,
758  {
759  typedef P1VTKFunction<GridView, Container> Function;
760  for (int c=0; c<ncomps; ++c) {
761  std::stringstream compName;
762  compName << name;
763  if (ncomps>1)
764  compName << "[" << c << "]";
765  VTKFunction* p = new Function(gridView_, v, compName.str(), ncomps, c, prec);
766  addVertexData(std::shared_ptr< const VTKFunction >(p));
767  }
768  }
769 
771  void clear ()
772  {
773  celldata.clear();
774  vertexdata.clear();
775  }
776 
779  { return coordPrec; }
780 
782  virtual ~VTKWriter ()
783  {
784  this->clear();
785  }
786 
798  std::string write ( const std::string &name,
799  VTK::OutputType type = VTK::ascii )
800  {
801  return write( name, type, gridView_.comm().rank(), gridView_.comm().size() );
802  }
803 
830  std::string pwrite ( const std::string & name, const std::string & path, const std::string & extendpath,
831  VTK::OutputType type = VTK::ascii )
832  {
833  return pwrite( name, path, extendpath, type, gridView_.comm().rank(), gridView_.comm().size() );
834  }
835 
836  protected:
838 
849  std::string getParallelPieceName(const std::string& name,
850  const std::string& path,
851  int commRank, int commSize) const
852  {
853  std::ostringstream s;
854  if(path.size() > 0) {
855  s << path;
856  if(path[path.size()-1] != '/')
857  s << '/';
858  }
859  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
860  s << 'p' << std::setw(4) << std::setfill('0') << commRank << '-';
861  s << name;
862  if(GridView::dimension > 1)
863  s << ".vtu";
864  else
865  s << ".vtp";
866  return s.str();
867  }
868 
870 
880  std::string getParallelHeaderName(const std::string& name,
881  const std::string& path,
882  int commSize) const
883  {
884  std::ostringstream s;
885  if(path.size() > 0) {
886  s << path;
887  if(path[path.size()-1] != '/')
888  s << '/';
889  }
890  s << 's' << std::setw(4) << std::setfill('0') << commSize << '-';
891  s << name;
892  if(GridView::dimension > 1)
893  s << ".pvtu";
894  else
895  s << ".pvtp";
896  return s.str();
897  }
898 
900 
912  std::string getSerialPieceName(const std::string& name,
913  const std::string& path) const
914  {
915  static const std::string extension =
916  GridView::dimension == 1 ? ".vtp" : ".vtu";
917 
918  return concatPaths(path, name+extension);
919  }
920 
936  std::string write ( const std::string &name,
937  VTK::OutputType type,
938  const int commRank,
939  const int commSize )
940  {
941  // in the parallel case, just use pwrite, it has all the necessary
942  // stuff, so we don't need to reimplement it here.
943  if(commSize > 1)
944  return pwrite(name, "", "", type, commRank, commSize);
945 
946  // make data mode visible to private functions
947  outputtype = type;
948 
949  // generate filename for process data
950  std::string pieceName = getSerialPieceName(name, "");
951 
952  // write process data
953  std::ofstream file;
954  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
955  std::ios_base::eofbit);
956  // check if file can be opened
957  try {
958  file.open( pieceName.c_str(), std::ios::binary );
959  }
960  catch(...) {
961  std::cerr << "Filename: " << pieceName << " could not be opened" << std::endl;
962  throw;
963  }
964  if (! file.is_open())
965  DUNE_THROW(IOError, "Could not write to piece file " << pieceName);
966  writeDataFile( file );
967  file.close();
968 
969  return pieceName;
970  }
971 
973 
996  std::string pwrite(const std::string& name, const std::string& path,
997  const std::string& extendpath,
998  VTK::OutputType ot, const int commRank,
999  const int commSize )
1000  {
1001  // make data mode visible to private functions
1002  outputtype=ot;
1003 
1004  // do some magic because paraview can only cope with relative paths to piece files
1005  std::ofstream file;
1006  file.exceptions(std::ios_base::badbit | std::ios_base::failbit |
1007  std::ios_base::eofbit);
1008  std::string piecepath = concatPaths(path, extendpath);
1009  std::string relpiecepath = relativePath(path, piecepath);
1010 
1011  // write this processes .vtu/.vtp piece file
1012  std::string fullname = getParallelPieceName(name, piecepath, commRank,
1013  commSize);
1014  // check if file can be opened
1015  try {
1016  file.open(fullname.c_str(),std::ios::binary);
1017  }
1018  catch(...) {
1019  std::cerr << "Filename: " << fullname << " could not be opened" << std::endl;
1020  throw;
1021  }
1022  if (! file.is_open())
1023  DUNE_THROW(IOError, "Could not write to piecefile file " << fullname);
1024  writeDataFile(file);
1025  file.close();
1026  gridView_.comm().barrier();
1027 
1028  // if we are rank 0, write .pvtu/.pvtp parallel header
1029  fullname = getParallelHeaderName(name, path, commSize);
1030  if( commRank ==0 )
1031  {
1032  file.open(fullname.c_str());
1033  if (! file.is_open())
1034  DUNE_THROW(IOError, "Could not write to parallel file " << fullname);
1035  writeParallelHeader(file,name,relpiecepath, commSize );
1036  file.close();
1037  }
1038  gridView_.comm().barrier();
1039  return fullname;
1040  }
1041 
1042  private:
1044 
1061  void writeParallelHeader(std::ostream& s, const std::string& piecename,
1062  const std::string& piecepath, const int commSize)
1063  {
1064  VTK::FileType fileType =
1065  (n == 1) ? VTK::polyData : VTK::unstructuredGrid;
1066 
1067  VTK::PVTUWriter writer(s, fileType);
1068 
1069  writer.beginMain();
1070 
1071  // PPointData
1072  {
1073  std::string scalars, vectors;
1074  std::tie(scalars,vectors) = getDataNames(vertexdata);
1075  writer.beginPointData(scalars, vectors);
1076  }
1077  for (auto it = vertexdata.begin(),
1078  end = vertexdata.end();
1079  it != end;
1080  ++it)
1081  {
1082  unsigned writecomps = it->fieldInfo().size();
1083  if(writecomps == 2) writecomps = 3;
1084  writer.addArray(it->name(), writecomps, it->fieldInfo().precision());
1085  }
1086  writer.endPointData();
1087 
1088  // PCellData
1089  {
1090  std::string scalars, vectors;
1091  std::tie(scalars,vectors) = getDataNames(celldata);
1092  writer.beginCellData(scalars, vectors);
1093  }
1094  for (auto it = celldata.begin(),
1095  end = celldata.end();
1096  it != end;
1097  ++it)
1098  {
1099  unsigned writecomps = it->fieldInfo().size();
1100  if(writecomps == 2) writecomps = 3;
1101  writer.addArray(it->name(), writecomps, it->fieldInfo().precision());
1102  }
1103  writer.endCellData();
1104 
1105  // PPoints
1106  writer.beginPoints();
1107  writer.addArray("Coordinates", 3, coordPrec);
1108  writer.endPoints();
1109 
1110  // Pieces
1111  for( int i = 0; i < commSize; ++i )
1112  {
1113  const std::string& fullname = getParallelPieceName(piecename,
1114  piecepath, i,
1115  commSize);
1116  writer.addPiece(fullname);
1117  }
1118 
1119  writer.endMain();
1120  }
1121 
1123  void writeDataFile (std::ostream& s)
1124  {
1125  VTK::FileType fileType =
1126  (n == 1) ? VTK::polyData : VTK::unstructuredGrid;
1127 
1128  VTK::VTUWriter writer(s, outputtype, fileType);
1129 
1130  // Grid characteristics
1131  vertexmapper = new VertexMapper( gridView_, mcmgVertexLayout() );
1132  if (datamode == VTK::conforming)
1133  {
1134  number.resize(vertexmapper->size());
1135  for (std::vector<int>::size_type i=0; i<number.size(); i++) number[i] = -1;
1136  }
1138 
1139  writer.beginMain(ncells, nvertices);
1140  writeAllData(writer);
1141  writer.endMain();
1142 
1143  // write appended binary data section
1144  if(writer.beginAppended())
1145  writeAllData(writer);
1146  writer.endAppended();
1147 
1148  delete vertexmapper; number.clear();
1149  }
1150 
1151  void writeAllData(VTK::VTUWriter& writer) {
1152  // PointData
1153  writeVertexData(writer);
1154 
1155  // CellData
1156  writeCellData(writer);
1157 
1158  // Points
1159  writeGridPoints(writer);
1160 
1161  // Cells
1162  writeGridCells(writer);
1163  }
1164 
1165  protected:
1166  std::string getFormatString() const
1167  {
1168  if (outputtype==VTK::ascii)
1169  return "ascii";
1170  if (outputtype==VTK::base64)
1171  return "binary";
1173  return "appended";
1175  return "appended";
1176  DUNE_THROW(IOError, "VTKWriter: unsupported OutputType" << outputtype);
1177  }
1178 
1179  std::string getTypeString() const
1180  {
1181  if (n==1)
1182  return "PolyData";
1183  else
1184  return "UnstructuredGrid";
1185  }
1186 
1188  virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
1189  {
1190  nvertices = 0;
1191  ncells = 0;
1192  ncorners = 0;
1193  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1194  {
1195  ncells++;
1196  // because of the use of vertexmapper->map(), this iteration must be
1197  // in the order of Dune's numbering.
1198  const int subEntities = it->subEntities(n);
1199  for (int i=0; i<subEntities; ++i)
1200  {
1201  ncorners++;
1202  if (datamode == VTK::conforming)
1203  {
1204  int alpha = vertexmapper->subIndex(*it,i,n);
1205  if (number[alpha]<0)
1206  number[alpha] = nvertices++;
1207  }
1208  else
1209  {
1210  nvertices++;
1211  }
1212  }
1213  }
1214  }
1215 
1216  template<typename T>
1217  std::tuple<std::string,std::string> getDataNames(const T& data) const
1218  {
1219  std::string scalars = "";
1220  for (auto it = data.begin(),
1221  end = data.end();
1222  it != end;
1223  ++it)
1224  if (it->fieldInfo().type() == VTK::FieldInfo::Type::scalar)
1225  {
1226  scalars = it->name();
1227  break;
1228  }
1229 
1230  std::string vectors = "";
1231  for (auto it = data.begin(),
1232  end = data.end();
1233  it != end;
1234  ++it)
1235  if (it->fieldInfo().type() == VTK::FieldInfo::Type::vector)
1236  {
1237  vectors = it->name();
1238  break;
1239  }
1240  return std::make_tuple(scalars,vectors);
1241  }
1242 
1243  template<typename Data, typename Iterator>
1244  void writeData(VTK::VTUWriter& writer, const Data& data, const Iterator begin, const Iterator end, int nentries)
1245  {
1246  for (auto it = data.begin(),
1247  iend = data.end();
1248  it != iend;
1249  ++it)
1250  {
1251  const auto& f = *it;
1252  VTK::FieldInfo fieldInfo = f.fieldInfo();
1253  std::size_t writecomps = fieldInfo.size();
1254  switch (fieldInfo.type())
1255  {
1257  break;
1259  // vtk file format: a vector data always should have 3 comps (with
1260  // 3rd comp = 0 in 2D case)
1261  if (writecomps > 3)
1262  DUNE_THROW(IOError,"Cannot write VTK vectors with more than 3 components (components was " << writecomps << ")");
1263  writecomps = 3;
1264  break;
1266  DUNE_THROW(NotImplemented,"VTK output for tensors not implemented yet");
1267  }
1268  std::shared_ptr<VTK::DataArrayWriter> p
1269  (writer.makeArrayWriter(f.name(), writecomps, nentries, fieldInfo.precision()));
1270  if(!p->writeIsNoop())
1271  for (Iterator eit = begin; eit!=end; ++eit)
1272  {
1273  const Entity & e = *eit;
1274  f.bind(e);
1275  f.write(eit.position(),*p);
1276  f.unbind();
1277  // vtk file format: a vector data always should have 3 comps
1278  // (with 3rd comp = 0 in 2D case)
1279  for (std::size_t j=fieldInfo.size(); j < writecomps; ++j)
1280  p->write(0.0);
1281  }
1282  }
1283  }
1284 
1286  virtual void writeCellData(VTK::VTUWriter& writer)
1287  {
1288  if(celldata.size() == 0)
1289  return;
1290 
1291  std::string scalars, vectors;
1292  std::tie(scalars,vectors) = getDataNames(celldata);
1293 
1294  writer.beginCellData(scalars, vectors);
1296  writer.endCellData();
1297  }
1298 
1300  virtual void writeVertexData(VTK::VTUWriter& writer)
1301  {
1302  if(vertexdata.size() == 0)
1303  return;
1304 
1305  std::string scalars, vectors;
1306  std::tie(scalars,vectors) = getDataNames(vertexdata);
1307 
1308  writer.beginPointData(scalars, vectors);
1310  writer.endPointData();
1311  }
1312 
1314  virtual void writeGridPoints(VTK::VTUWriter& writer)
1315  {
1316  writer.beginPoints();
1317 
1318  std::shared_ptr<VTK::DataArrayWriter> p
1319  (writer.makeArrayWriter("Coordinates", 3, nvertices, coordPrec));
1320  if(!p->writeIsNoop()) {
1321  VertexIterator vEnd = vertexEnd();
1322  for (VertexIterator vit=vertexBegin(); vit!=vEnd; ++vit)
1323  {
1324  int dimw=w;
1325  for (int j=0; j<std::min(dimw,3); j++)
1326  p->write((*vit).geometry().corner(vit.localindex())[j]);
1327  for (int j=std::min(dimw,3); j<3; j++)
1328  p->write(0.0);
1329  }
1330  }
1331  // free the VTK::DataArrayWriter before touching the stream
1332  p.reset();
1333 
1334  writer.endPoints();
1335  }
1336 
1338  virtual void writeGridCells(VTK::VTUWriter& writer)
1339  {
1340  writer.beginCells();
1341 
1342  // connectivity
1343  {
1344  std::shared_ptr<VTK::DataArrayWriter> p1
1345  (writer.makeArrayWriter("connectivity", 1, ncorners, VTK::Precision::int32));
1346  if(!p1->writeIsNoop())
1347  for (CornerIterator it=cornerBegin(); it!=cornerEnd(); ++it)
1348  p1->write(it.id());
1349  }
1350 
1351  // offsets
1352  {
1353  std::shared_ptr<VTK::DataArrayWriter> p2
1354  (writer.makeArrayWriter("offsets", 1, ncells, VTK::Precision::int32));
1355  if(!p2->writeIsNoop()) {
1356  int offset = 0;
1357  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1358  {
1359  offset += it->subEntities(n);
1360  p2->write(offset);
1361  }
1362  }
1363  }
1364 
1365  // types
1366  if (n>1)
1367  {
1368  {
1369  std::shared_ptr<VTK::DataArrayWriter> p3
1370  (writer.makeArrayWriter("types", 1, ncells, VTK::Precision::uint8));
1371 
1372  if(!p3->writeIsNoop())
1373  {
1374  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
1375  {
1376  int vtktype = VTK::geometryType(it->type());
1377  p3->write(vtktype);
1378  }
1379  }
1380  }
1381 
1382 
1383  // if polyhedron cells found also cell faces need to be written
1384  if( polyhedralCellsPresent_ )
1385  {
1386  writeCellFaces( writer );
1387  }
1388  }
1389 
1390  writer.endCells();
1391  }
1392 
1393  protected:
1395  {
1396  // check if polyhedron cells are present
1397  for( const auto& geomType : gridView_.indexSet().types( 0 ) )
1398  {
1399  if( VTK::geometryType( geomType ) == VTK::polyhedron )
1400  {
1401  return true;
1402  }
1403  }
1404  return false;
1405  }
1406 
1408  virtual void writeCellFaces(VTK::VTUWriter& writer)
1409  {
1410  if( ! faceVertices_ )
1411  {
1412  faceVertices_.reset( new std::pair< std::vector<int>, std::vector<int> > () );
1413  // fill face vertex structure
1415  faceVertices_->first, faceVertices_->second );
1416  }
1417 
1418  std::vector< int >& faces = faceVertices_->first;
1419  std::vector< int >& faceOffsets = faceVertices_->second;
1420  assert( int(faceOffsets.size()) == ncells );
1421 
1422  {
1423  std::shared_ptr<VTK::DataArrayWriter> p4
1424  (writer.makeArrayWriter("faces", 1, faces.size(), VTK::Precision::int32));
1425  if(!p4->writeIsNoop())
1426  {
1427  for( const auto& face : faces )
1428  p4->write( face );
1429  }
1430  }
1431 
1432  {
1433  std::shared_ptr<VTK::DataArrayWriter> p5
1434  (writer.makeArrayWriter("faceoffsets", 1, ncells, VTK::Precision::int32));
1435  if(!p5->writeIsNoop())
1436  {
1437  for( const auto& offset : faceOffsets )
1438  p5->write( offset );
1439 
1440  // clear face vertex structure
1441  faceVertices_.reset();
1442  }
1443  }
1444  }
1445 
1446  template <class CornerIterator, class IndexSet, class T>
1448  const CornerIterator end,
1449  const IndexSet& indexSet,
1450  std::vector<T>& faces,
1451  std::vector<T>& faceOffsets )
1452  {
1453  if( n == 3 && it != end )
1454  {
1455  // clear output arrays
1456  faces.clear();
1457  faces.reserve( 15 * ncells );
1458  faceOffsets.clear();
1459  faceOffsets.reserve( ncells );
1460 
1461  int offset = 0;
1462 
1463  Cell element = *it;
1464  int elIndex = indexSet.index( element );
1465  std::vector< T > vertices;
1466  vertices.reserve( 30 );
1467  for( ; it != end; ++it )
1468  {
1469  const Cell& cell = *it ;
1470  const int cellIndex = indexSet.index( cell ) ;
1471  if( elIndex != cellIndex )
1472  {
1473  fillFacesForElement( element, indexSet, vertices, offset, faces, faceOffsets );
1474 
1475  vertices.clear();
1476  element = cell ;
1477  elIndex = cellIndex ;
1478  }
1479  vertices.push_back( it.id() );
1480  }
1481 
1482  // fill faces for last element
1483  fillFacesForElement( element, indexSet, vertices, offset, faces, faceOffsets );
1484  }
1485  }
1486 
1487  template <class Entity, class IndexSet, class T>
1488  static void fillFacesForElement( const Entity& element,
1489  const IndexSet& indexSet,
1490  const std::vector<T>& vertices,
1491  T& offset,
1492  std::vector<T>& faces,
1493  std::vector<T>& faceOffsets )
1494  {
1495  const int dim = n;
1496 
1497  std::map< T, T > vxMap;
1498 
1499  // get number of local faces
1500  const int nVertices = element.subEntities( dim );
1501  for( int vx = 0; vx < nVertices; ++ vx )
1502  {
1503  const int vxIdx = indexSet.subIndex( element, vx, dim );
1504  vxMap[ vxIdx ] = vertices[ vx ];
1505  }
1506 
1507  // get number of local faces
1508  const int nFaces = element.subEntities( 1 );
1509  // store number of faces for current element
1510  faces.push_back( nFaces );
1511  ++offset;
1512  // extract each face as a set of vertex indices
1513  for( int fce = 0; fce < nFaces; ++ fce )
1514  {
1515  // obtain face
1516  const auto face = element.template subEntity< 1 > ( fce );
1517 
1518  // get all vertex indices from current face
1519  const int nVxFace = face.subEntities( dim );
1520  faces.push_back( nVxFace );
1521  ++offset ;
1522  for( int i=0; i<nVxFace; ++i )
1523  {
1524  const T vxIndex = indexSet.subIndex( face, i, dim );
1525  assert( vxMap.find( vxIndex ) != vxMap.end() );
1526  faces.push_back( vxMap[ vxIndex ] );
1527  ++offset ;
1528  }
1529  }
1530 
1531  // store face offset for each element
1532  faceOffsets.push_back( offset );
1533  }
1534 
1535  protected:
1536  // the list of registered functions
1537  std::list<VTKLocalFunction> celldata;
1538  std::list<VTKLocalFunction> vertexdata;
1539 
1540  // the grid
1542 
1543  // temporary grid information
1544  int ncells;
1547  private:
1548  VertexMapper* vertexmapper;
1549  // in conforming mode, for each vertex id (as obtained by vertexmapper)
1550  // hold its number in the iteration order (VertexIterator)
1551  std::vector<int> number;
1552  VTK::DataMode datamode;
1553  VTK::Precision coordPrec;
1554 
1555  // true if polyhedral cells are present in the grid
1556  const bool polyhedralCellsPresent_;
1557 
1558  // pointer holding face vertex connectivity if needed
1559  std::shared_ptr< std::pair< std::vector<int>, std::vector<int> > > faceVertices_;
1560 
1561  protected:
1563  };
1564 
1565 }
1566 
1567 #endif
Mapper for multiple codim and multiple geometry types.
Common stuff for the VTKWriter.
Data array writers for the VTKWriter.
Functions for VTK output.
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
Traits ::IndexSet IndexSet
type of the index set
Definition: common/gridview.hh:80
Grid::ctype ctype
type used for coordinates in grid
Definition: common/gridview.hh:124
const CollectiveCommunication & comm() const
obtain collective communication object
Definition: common/gridview.hh:246
Traits ::Grid Grid
type of the grid
Definition: common/gridview.hh:77
const IndexSet & indexSet() const
obtain the index set
Definition: common/gridview.hh:172
@ dimensionworld
The dimension of the world the grid lives in.
Definition: common/gridview.hh:131
@ dimension
The dimension of the grid.
Definition: common/gridview.hh:127
MCMGLayout mcmgVertexLayout()
layout for vertices (dim-0 entities)
Definition: mcmgmapper.hh:160
Include standard header files.
Definition: agrid.hh:59
int min(const DofVectorPointer< int > &dofVector)
Definition: dofvector.hh:346
Precision
which precision to use when writing out data to vtk files
Definition: common.hh:319
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:42
@ ascii
Output to the file is in ascii.
Definition: common.hh:44
@ appendedraw
Output is to the file is appended raw binary.
Definition: common.hh:48
@ appendedbase64
Output is to the file is appended base64 binary.
Definition: common.hh:50
@ base64
Output to the file is inline base64 binary.
Definition: common.hh:46
int renumber(const Dune::GeometryType &t, int i)
renumber VTK <-> Dune
Definition: common.hh:234
FileType
which type of VTK file to write
Definition: common.hh:300
@ polyData
for .vtp files (PolyData)
Definition: common.hh:302
@ unstructuredGrid
for .vtu files (UnstructuredGrid)
Definition: common.hh:304
DataMode
Whether to produce conforming or non-conforming output.
Definition: common.hh:66
@ conforming
Output conforming data.
Definition: common.hh:72
@ nonconforming
Output non-conforming data.
Definition: common.hh:80
GeometryType geometryType(const Dune::GeometryType &t)
mapping from GeometryType to VTKGeometryType
Definition: common.hh:199
@ polyhedron
Definition: common.hh:190
Grid view abstract base class.
Definition: common/gridview.hh:60
Implementation class for a multiple codim and multiple geometry type mapper.
Definition: mcmgmapper.hh:200
Index subIndex(const typename GV::template Codim< 0 >::Entity &e, int i, unsigned int codim) const
Map subentity of codim 0 entity to starting index in array for dof block.
Definition: mcmgmapper.hh:268
size_type size() const
Return total number of entities in the entity set managed by the mapper.
Definition: mcmgmapper.hh:287
Descriptor struct for VTK fields.
Definition: common.hh:376
std::size_t size() const
The number of components in the data field.
Definition: common.hh:412
Precision precision() const
The precision used for the output of the data field.
Definition: common.hh:418
@ tensor
tensor field (always 3x3)
@ vector
vector-valued field (always 3D, will be padded if necessary)
Type type() const
The type of the data field.
Definition: common.hh:406
std::string name() const
The name of the data field.
Definition: common.hh:400
base class for data array writers
Definition: dataarraywriter.hh:54
void write(T data)
write one element of data
Definition: dataarraywriter.hh:67
A base class for grid functions with any return type and dimension.
Definition: function.hh:40
Take a vector and interpret it as cell data for the VTKWriter.
Definition: function.hh:95
Take a vector and interpret it as point data for the VTKWriter.
Definition: function.hh:203
Dump a .vtu/.vtp files contents to a stream.
Definition: pvtuwriter.hh:60
Writer for the ouput of grid functions in the vtk format.
Definition: vtksequencewriter.hh:28
Base class to write pvd-files which contains a list of all collected vtk-files.
Definition: vtksequencewriterbase.hh:32
Writer for the ouput of grid functions in the vtk format.
Definition: vtkwriter.hh:91
void addCellData(const Container &v, const std::string &name, int ncomps=1, VTK::Precision prec=VTK::Precision::float32)
Add a grid function (represented by container) that lives on the cells of the grid to the visualizati...
Definition: vtkwriter.hh:691
CornerIterator cornerEnd() const
Definition: vtkwriter.hh:616
void clear()
clear list of registered functions
Definition: vtkwriter.hh:771
std::string write(const std::string &name, VTK::OutputType type=VTK::ascii)
write output (interface might change later)
Definition: vtkwriter.hh:798
VertexIterator vertexBegin() const
Definition: vtkwriter.hh:504
std::string getTypeString() const
Definition: vtkwriter.hh:1179
std::string getParallelHeaderName(const std::string &name, const std::string &path, int commSize) const
return name of a parallel header file
Definition: vtkwriter.hh:880
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:709
Dune::VTKFunction< GridView > VTKFunction
Definition: vtkwriter.hh:143
CellIterator cellEnd() const
Definition: vtkwriter.hh:398
std::list< VTKLocalFunction > vertexdata
Definition: vtkwriter.hh:1538
std::tuple< std::string, std::string > getDataNames(const T &data) const
Definition: vtkwriter.hh:1217
CornerIterator cornerBegin() const
Definition: vtkwriter.hh:609
std::string getSerialPieceName(const std::string &name, const std::string &path) const
return name of a serial piece file
Definition: vtkwriter.hh:912
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:645
std::string getFormatString() const
Definition: vtkwriter.hh:1166
bool checkForPolyhedralCells() const
Definition: vtkwriter.hh:1394
void addVertexData(F &&f, VTK::FieldInfo vtkFieldInfo)
Add a function by sampling it on the grid vertices.
Definition: vtkwriter.hh:734
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: vtkwriter.hh:1286
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:849
CellIterator cellBegin() const
Definition: vtkwriter.hh:393
virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
count the vertices, cells and corners
Definition: vtkwriter.hh:1188
VTK::OutputType outputtype
Definition: vtkwriter.hh:1562
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: vtkwriter.hh:1338
GridView gridView_
Definition: vtkwriter.hh:1541
virtual void writeCellFaces(VTK::VTUWriter &writer)
write the connectivity array
Definition: vtkwriter.hh:1408
void fillFaceVertices(CornerIterator it, const CornerIterator end, const IndexSet &indexSet, std::vector< T > &faces, std::vector< T > &faceOffsets)
Definition: vtkwriter.hh:1447
std::list< VTKLocalFunction > celldata
Definition: vtkwriter.hh:1537
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:936
VTK::Precision coordPrecision() const
get the precision with which coordinates are written out
Definition: vtkwriter.hh:778
std::list< VTKLocalFunction >::const_iterator FunctionIterator
Definition: vtkwriter.hh:372
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: vtkwriter.hh:1314
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: vtkwriter.hh:1300
int nvertices
Definition: vtkwriter.hh:1545
void addCellData(F &&f, VTK::FieldInfo vtkFieldInfo)
Add a function by sampling it on the element centers.
Definition: vtkwriter.hh:670
void addVertexData(const Container &v, const std::string &name, int ncomps=1, VTK::Precision prec=VTK::Precision::float32)
Add a grid function (represented by container) that lives on the vertices of the grid to the visualiz...
Definition: vtkwriter.hh:756
virtual ~VTKWriter()
destructor
Definition: vtkwriter.hh:782
static void fillFacesForElement(const Entity &element, const IndexSet &indexSet, const std::vector< T > &vertices, T &offset, std::vector< T > &faces, std::vector< T > &faceOffsets)
Definition: vtkwriter.hh:1488
void writeData(VTK::VTUWriter &writer, const Data &data, const Iterator begin, const Iterator end, int nentries)
Definition: vtkwriter.hh:1244
int ncells
Definition: vtkwriter.hh:1544
VertexIterator vertexEnd() const
Definition: vtkwriter.hh:511
VTKWriter(const GridView &gridView, VTK::DataMode dm=VTK::conforming, VTK::Precision coordPrecision=VTK::Precision::float32)
Construct a VTKWriter working on a specific GridView.
Definition: vtkwriter.hh:632
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:996
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:830
int ncorners
Definition: vtkwriter.hh:1546
Type erasure wrapper for VTK data sets.
Definition: vtkwriter.hh:152
void unbind() const
Unbind the data set from the currently bound entity.
Definition: vtkwriter.hh:356
VTKLocalFunction(F &&f, VTK::FieldInfo fieldInfo)
Construct a VTKLocalFunction for a dune-functions style LocalFunction.
Definition: vtkwriter.hh:303
std::string name() const
Returns the name of the data set.
Definition: vtkwriter.hh:338
const VTK::FieldInfo & fieldInfo() const
Returns the VTK::FieldInfo for the data set.
Definition: vtkwriter.hh:344
VTK::FieldInfo _fieldInfo
Definition: vtkwriter.hh:368
VTK::DataArrayWriter Writer
Definition: vtkwriter.hh:156
void bind(const Entity &e) const
Bind the data set to grid entity e.
Definition: vtkwriter.hh:350
VTKLocalFunction(const std::shared_ptr< const VTKFunction > &vtkFunctionPtr)
Construct a VTKLocalFunction for a legacy VTKFunction.
Definition: vtkwriter.hh:327
std::shared_ptr< FunctionWrapperBase > _f
Definition: vtkwriter.hh:367
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:362
Base class for polymorphic container of underlying data set.
Definition: vtkwriter.hh:160
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 ~FunctionWrapperBase()
Definition: vtkwriter.hh:174
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:183
typename std::decay< F >::type Function
Definition: vtkwriter.hh:184
FunctionWrapper(F_ &&f)
Definition: vtkwriter.hh:187
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:201
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:196
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:191
Type erasure implementation for C++ functions, i.e., functions that can be evaluated in global coordi...
Definition: vtkwriter.hh:231
GlobalFunctionWrapper(F_ &&f)
Definition: vtkwriter.hh:235
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:245
typename std::decay< F >::type Function
Definition: vtkwriter.hh:232
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:250
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:240
Type erasure implementation for legacy VTKFunctions.
Definition: vtkwriter.hh:272
virtual void unbind()
Unbind data set from current grid entity - mostly here for performance and symmetry reasons.
Definition: vtkwriter.hh:283
VTKFunctionWrapper(const std::shared_ptr< const VTKFunction > &f)
Definition: vtkwriter.hh:273
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:288
virtual void bind(const Entity &e)
Bind data set to grid entity - must be called before evaluating (i.e. calling write())
Definition: vtkwriter.hh:278
Iterator over the grids elements.
Definition: vtkwriter.hh:381
const FieldVector< DT, n > position() const
Definition: vtkwriter.hh:387
CellIterator(const GridCellIterator &x)
construct a CellIterator from the gridview's Iterator.
Definition: vtkwriter.hh:384
Iterate over the grid's vertices.
Definition: vtkwriter.hh:420
FieldVector< DT, n > position() const
position of vertex inside the entity
Definition: vtkwriter.hh:497
VertexIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm)
Definition: vtkwriter.hh:453
void basicIncrement()
Definition: vtkwriter.hh:436
void increment()
Definition: vtkwriter.hh:464
EntityReference dereference() const
Definition: vtkwriter.hh:487
bool equals(const VertexIterator &cit) const
Definition: vtkwriter.hh:481
int localindex() const
index of vertex within the entity, in Dune-numbering
Definition: vtkwriter.hh:492
Iterate over the elements' corners.
Definition: vtkwriter.hh:535
void increment()
Definition: vtkwriter.hh:562
CornerIterator(const GridCellIterator &x, const GridCellIterator &end, const VTK::DataMode &dm, const VertexMapper &vm, const std::vector< int > &num)
Definition: vtkwriter.hh:554
int id() const
Process-local consecutive zero-starting vertex id.
Definition: vtkwriter.hh:593
EntityReference dereference() const
Definition: vtkwriter.hh:584
bool equals(const CornerIterator &cit) const
Definition: vtkwriter.hh:578
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
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
DataArrayWriter * makeArrayWriter(const std::string &name, unsigned ncomps, unsigned nitems, Precision prec)
acquire a DataArrayWriter
Definition: vtuwriter.hh:378
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