dune-grid  2.8.0
subsamplingvtkwriter.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_SUBSAMPLINGVTKWRITER_HH
5 #define DUNE_SUBSAMPLINGVTKWRITER_HH
6 
7 #include <ostream>
8 #include <memory>
9 
10 #include <dune/common/indent.hh>
11 #include <dune/geometry/type.hh>
12 #include <dune/geometry/virtualrefinement.hh>
15 
22 namespace Dune
23 {
35  template< class GridView >
37  : public VTKWriter<GridView>
38  {
39  typedef VTKWriter<GridView> Base;
40  enum { dim = GridView::dimension };
41  enum { dimw = GridView::dimensionworld };
42  typedef typename GridView::Grid::ctype ctype;
43  typedef typename GridView::template Codim< 0 >::Entity Entity;
44  typedef VirtualRefinement<dim, ctype> Refinement;
45  typedef typename Refinement::IndexVector IndexVector;
46  typedef typename Refinement::ElementIterator SubElementIterator;
47  typedef typename Refinement::VertexIterator SubVertexIterator;
48 
49  typedef typename Base::CellIterator CellIterator;
50  typedef typename Base::FunctionIterator FunctionIterator;
51  using Base::cellBegin;
52  using Base::cellEnd;
53  using Base::celldata;
54  using Base::ncells;
55  using Base::ncorners;
56  using Base::nvertices;
57  using Base::outputtype;
58  using Base::vertexBegin;
59  using Base::vertexEnd;
60  using Base::vertexdata;
61 
62  public:
78  explicit SubsamplingVTKWriter (const GridView &gridView,
79  Dune::RefinementIntervals intervals_, bool coerceToSimplex_ = false,
81  : Base(gridView, VTK::nonconforming, coordPrecision)
82  , intervals(intervals_), coerceToSimplex(coerceToSimplex_)
83  {
84  if(intervals_.intervals() < 1) {
85  DUNE_THROW(Dune::IOError,"SubsamplingVTKWriter: Refinement intervals must be larger than zero! (One interval means no subsampling)");
86  }
87  }
88 
89  private:
90  GeometryType subsampledGeometryType(GeometryType geometryType)
91  {
92  return (geometryType.isCube() && !coerceToSimplex ? geometryType : GeometryTypes::simplex(dim));
93  }
94 
95  template<typename SubIterator>
96  struct IteratorSelector
97  {};
98 
99  SubElementIterator refinementBegin(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubElementIterator>)
100  {
101  return refinement.eBegin(intervals);
102  }
103 
104  SubVertexIterator refinementBegin(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubVertexIterator>)
105  {
106  return refinement.vBegin(intervals);
107  }
108 
109  SubElementIterator refinementEnd(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubElementIterator>)
110  {
111  return refinement.eEnd(intervals);
112  }
113 
114  SubVertexIterator refinementEnd(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubVertexIterator>)
115  {
116  return refinement.vEnd(intervals);
117  }
118 
119  template<typename Data, typename Iterator, typename SubIterator>
120  void writeData(VTK::VTUWriter& writer, const Data& data, const Iterator begin, const Iterator end, int nentries, IteratorSelector<SubIterator> sis)
121  {
122  for (auto it = data.begin(),
123  iend = data.end();
124  it != iend;
125  ++it)
126  {
127  const auto& f = *it;
128  VTK::FieldInfo fieldInfo = f.fieldInfo();
129  std::size_t writecomps = fieldInfo.size();
130  switch (fieldInfo.type())
131  {
133  break;
135  // vtk file format: a vector data always should have 3 comps (with
136  // 3rd comp = 0 in 2D case)
137  if (writecomps > 3)
138  DUNE_THROW(IOError,"Cannot write VTK vectors with more than 3 components (components was " << writecomps << ")");
139  writecomps = 3;
140  break;
142  DUNE_THROW(NotImplemented,"VTK output for tensors not implemented yet");
143  }
144  std::shared_ptr<VTK::DataArrayWriter> p
145  (writer.makeArrayWriter(f.name(), writecomps, nentries, fieldInfo.precision()));
146  if(!p->writeIsNoop())
147  for (Iterator eit = begin; eit!=end; ++eit)
148  {
149  const Entity & e = *eit;
150  f.bind(e);
151  Refinement &refinement =
152  buildRefinement<dim, ctype>(eit->type(),
153  subsampledGeometryType(eit->type()));
154  for(SubIterator sit = refinementBegin(refinement,intervals,sis),
155  send = refinementEnd(refinement,intervals,sis);
156  sit != send;
157  ++sit)
158  {
159  f.write(sit.coords(),*p);
160  // expand 2D-Vectors to 3D for VTK format
161  for(unsigned j = f.fieldInfo().size(); j < writecomps; j++)
162  p->write(0.0);
163  }
164  f.unbind();
165  }
166  }
167  }
168 
169 
170  protected:
172  virtual void countEntities(int &nvertices_, int &ncells_, int &ncorners_);
173 
175  virtual void writeCellData(VTK::VTUWriter& writer);
176 
178  virtual void writeVertexData(VTK::VTUWriter& writer);
179 
181  virtual void writeGridPoints(VTK::VTUWriter& writer);
182 
184  virtual void writeGridCells(VTK::VTUWriter& writer);
185 
186  public:
187  using Base::addVertexData;
188  using Base::addCellData;
189 
190  private:
191  // hide addVertexData -- adding raw data directly without a VTKFunction
192  // currently does not make sense for subsampled meshes, as the higher order
193  // information is missing. See FS#676.
194  template<class V>
195  void addVertexData (const V& v, const std::string &name, int ncomps=1);
196  template<class V>
197  void addCellData (const V& v, const std::string &name, int ncomps=1);
198 
199  Dune::RefinementIntervals intervals;
200  bool coerceToSimplex;
201  };
202 
204  template <class GridView>
205  void SubsamplingVTKWriter<GridView>::countEntities(int &nvertices_, int &ncells_, int &ncorners_)
206  {
207  nvertices_ = 0;
208  ncells_ = 0;
209  ncorners_ = 0;
210  for (CellIterator it=this->cellBegin(); it!=cellEnd(); ++it)
211  {
212  Refinement &refinement = buildRefinement<dim, ctype>(it->type(), subsampledGeometryType(it->type()));
213 
214  ncells_ += refinement.nElements(intervals);
215  nvertices_ += refinement.nVertices(intervals);
216  ncorners_ += refinement.nElements(intervals) * refinement.eBegin(intervals).vertexIndices().size();
217  }
218  }
219 
220 
222  template <class GridView>
224  {
225  if(celldata.size() == 0)
226  return;
227 
228  // Find the names of the first scalar and vector data fields.
229  // These will be marked as the default fields (the ones that ParaView shows
230  // when the file has just been opened).
231  std::string defaultScalarField, defaultVectorField;
232  std::tie(defaultScalarField, defaultVectorField) = this->getDataNames(celldata);
233 
234  writer.beginCellData(defaultScalarField, defaultVectorField);
235  writeData(writer,celldata,cellBegin(),cellEnd(),ncells,IteratorSelector<SubElementIterator>());
236  writer.endCellData();
237  }
238 
240  template <class GridView>
242  {
243  if(vertexdata.size() == 0)
244  return;
245 
246  // Find the names of the first scalar and vector data fields.
247  // These will be marked as the default fields (the ones that ParaView shows
248  // when the file has just been opened).
249  std::string defaultScalarField, defaultVectorField;
250  std::tie(defaultScalarField, defaultVectorField) = this->getDataNames(vertexdata);
251 
252  writer.beginPointData(defaultScalarField, defaultVectorField);
253  writeData(writer,vertexdata,cellBegin(),cellEnd(),nvertices,IteratorSelector<SubVertexIterator>());
254  writer.endPointData();
255  }
256 
258  template <class GridView>
260  {
261  writer.beginPoints();
262 
263  std::shared_ptr<VTK::DataArrayWriter> p
264  (writer.makeArrayWriter("Coordinates", 3, nvertices, this->coordPrecision()));
265  if(!p->writeIsNoop())
266  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
267  {
268  Refinement &refinement =
269  buildRefinement<dim, ctype>(i->type(),
270  subsampledGeometryType(i->type()));
271  for(SubVertexIterator sit = refinement.vBegin(intervals),
272  send = refinement.vEnd(intervals);
273  sit != send; ++sit)
274  {
275  FieldVector<ctype, dimw> coords = i->geometry().global(sit.coords());
276  for (int j=0; j<std::min(int(dimw),3); j++)
277  p->write(coords[j]);
278  for (int j=std::min(int(dimw),3); j<3; j++)
279  p->write(0.0);
280  }
281  }
282  // free the VTK::DataArrayWriter before touching the stream
283  p.reset();
284 
285  writer.endPoints();
286  }
287 
289  template <class GridView>
291  {
292  writer.beginCells();
293 
294  // connectivity
295  {
296  std::shared_ptr<VTK::DataArrayWriter> p1
297  (writer.makeArrayWriter("connectivity", 1, ncorners, VTK::Precision::int32));
298  // The offset within the index numbering
299  if(!p1->writeIsNoop()) {
300  int offset = 0;
301  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
302  {
303  GeometryType coercedToType = subsampledGeometryType(i->type());
304  Refinement &refinement =
305  buildRefinement<dim, ctype>(i->type(), coercedToType);
306  for(SubElementIterator sit = refinement.eBegin(intervals),
307  send = refinement.eEnd(intervals);
308  sit != send; ++sit)
309  {
310  IndexVector indices = sit.vertexIndices();
311  for(unsigned int ii = 0; ii < indices.size(); ++ii)
312  p1->write(offset+indices[VTK::renumber(coercedToType, ii)]);
313  }
314  offset += refinement.nVertices(intervals);
315  }
316  }
317  }
318 
319  // offsets
320  {
321  std::shared_ptr<VTK::DataArrayWriter> p2
322  (writer.makeArrayWriter("offsets", 1, ncells, VTK::Precision::int32));
323  if(!p2->writeIsNoop()) {
324  // The offset into the connectivity array
325  int offset = 0;
326  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
327  {
328  Refinement &refinement =
329  buildRefinement<dim, ctype>(i->type(),
330  subsampledGeometryType(i->type()));
331  unsigned int verticesPerCell =
332  refinement.eBegin(intervals).vertexIndices().size();
333  for(int element = 0; element < refinement.nElements(intervals);
334  ++element)
335  {
336  offset += verticesPerCell;
337  p2->write(offset);
338  }
339  }
340  }
341  }
342 
343  // types
344  if (dim>1)
345  {
346  std::shared_ptr<VTK::DataArrayWriter> p3
347  (writer.makeArrayWriter("types", 1, ncells, VTK::Precision::uint8));
348  if(!p3->writeIsNoop())
349  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
350  {
351  GeometryType coerceTo = subsampledGeometryType(it->type());
352  Refinement &refinement =
353  buildRefinement<dim, ctype>(it->type(), coerceTo);
354  int vtktype = VTK::geometryType(coerceTo);
355  for(int i = 0; i < refinement.nElements(intervals); ++i)
356  p3->write(vtktype);
357  }
358  }
359 
360  writer.endCells();
361  }
362 }
363 
364 #endif // DUNE_SUBSAMPLINGVTKWRITER_HH
Provides file i/o for the visualization toolkit.
@ dimensionworld
The dimension of the world the grid lives in.
Definition: common/gridview.hh:134
@ dimension
The dimension of the grid.
Definition: common/gridview.hh:130
Include standard header files.
Definition: agrid.hh:58
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:269
int renumber(const Dune::GeometryType &t, int i)
renumber VTK <-> Dune
Definition: common.hh:184
@ nonconforming
Output non-conforming data.
Definition: common.hh:79
GeometryType geometryType(const Dune::GeometryType &t)
mapping from GeometryType to VTKGeometryType
Definition: common.hh:149
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:130
Grid view abstract base class.
Definition: common/gridview.hh:63
@ tensor
tensor field (always 3x3)
@ vector
vector-valued field (always 3D, will be padded if necessary)
Writer for the output of subsampled grid functions in the vtk format.
Definition: subsamplingvtkwriter.hh:38
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:711
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:647
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: subsamplingvtkwriter.hh:259
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: subsamplingvtkwriter.hh:241
SubsamplingVTKWriter(const GridView &gridView, Dune::RefinementIntervals intervals_, bool coerceToSimplex_=false, VTK::Precision coordPrecision=VTK::Precision::float32)
Construct a SubsamplingVTKWriter working on a specific GridView.
Definition: subsamplingvtkwriter.hh:78
virtual void countEntities(int &nvertices_, int &ncells_, int &ncorners_)
count the vertices, cells and corners
Definition: subsamplingvtkwriter.hh:205
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: subsamplingvtkwriter.hh:223
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: subsamplingvtkwriter.hh:290
Writer for the ouput of grid functions in the vtk format.
Definition: vtkwriter.hh:93
VertexIterator vertexBegin() const
Definition: vtkwriter.hh:506
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:711
CellIterator cellEnd() const
Definition: vtkwriter.hh:400
std::list< VTKLocalFunction > vertexdata
Definition: vtkwriter.hh:1583
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:647
CellIterator cellBegin() const
Definition: vtkwriter.hh:395
VTK::OutputType outputtype
Definition: vtkwriter.hh:1607
std::list< VTKLocalFunction > celldata
Definition: vtkwriter.hh:1582
VTK::Precision coordPrecision() const
get the precision with which coordinates are written out
Definition: vtkwriter.hh:780
std::list< VTKLocalFunction >::const_iterator FunctionIterator
Definition: vtkwriter.hh:374
int nvertices
Definition: vtkwriter.hh:1590
int ncells
Definition: vtkwriter.hh:1589
VertexIterator vertexEnd() const
Definition: vtkwriter.hh:513
int ncorners
Definition: vtkwriter.hh:1591
Iterator over the grids elements.
Definition: vtkwriter.hh:383
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