dune-grid  2.7.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: Negative Subsampling " << intervals_.intervals() << " must not be used!");
86  }
87  }
100  DUNE_DEPRECATED_MSG("SubsampligVTKWriter(GV,int,bool) is deprecated, use SubsamplingVTKWriter(GV,Dune::refinement{Intervals|Levels}(int),bool)")
101  explicit SubsamplingVTKWriter (const GridView &gridView,
102  int level_, bool coerceToSimplex_ = false)
103  : SubsamplingVTKWriter(gridView, Dune::refinementIntervals(1<<level_), coerceToSimplex_)
104  { }
105 
106  private:
107  GeometryType subsampledGeometryType(GeometryType geometryType)
108  {
109  return (geometryType.isCube() && !coerceToSimplex ? geometryType : GeometryTypes::simplex(dim));
110  }
111 
112  template<typename SubIterator>
113  struct IteratorSelector
114  {};
115 
116  SubElementIterator refinementBegin(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubElementIterator>)
117  {
118  return refinement.eBegin(intervals);
119  }
120 
121  SubVertexIterator refinementBegin(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubVertexIterator>)
122  {
123  return refinement.vBegin(intervals);
124  }
125 
126  SubElementIterator refinementEnd(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubElementIterator>)
127  {
128  return refinement.eEnd(intervals);
129  }
130 
131  SubVertexIterator refinementEnd(const Refinement& refinement, Dune::RefinementIntervals intervals, IteratorSelector<SubVertexIterator>)
132  {
133  return refinement.vEnd(intervals);
134  }
135 
136  template<typename Data, typename Iterator, typename SubIterator>
137  void writeData(VTK::VTUWriter& writer, const Data& data, const Iterator begin, const Iterator end, int nentries, IteratorSelector<SubIterator> sis)
138  {
139  for (auto it = data.begin(),
140  iend = data.end();
141  it != iend;
142  ++it)
143  {
144  const auto& f = *it;
145  VTK::FieldInfo fieldInfo = f.fieldInfo();
146  std::size_t writecomps = fieldInfo.size();
147  switch (fieldInfo.type())
148  {
150  break;
152  // vtk file format: a vector data always should have 3 comps (with
153  // 3rd comp = 0 in 2D case)
154  if (writecomps > 3)
155  DUNE_THROW(IOError,"Cannot write VTK vectors with more than 3 components (components was " << writecomps << ")");
156  writecomps = 3;
157  break;
159  DUNE_THROW(NotImplemented,"VTK output for tensors not implemented yet");
160  }
161  std::shared_ptr<VTK::DataArrayWriter> p
162  (writer.makeArrayWriter(f.name(), writecomps, nentries, fieldInfo.precision()));
163  if(!p->writeIsNoop())
164  for (Iterator eit = begin; eit!=end; ++eit)
165  {
166  const Entity & e = *eit;
167  f.bind(e);
168  Refinement &refinement =
169  buildRefinement<dim, ctype>(eit->type(),
170  subsampledGeometryType(eit->type()));
171  for(SubIterator sit = refinementBegin(refinement,intervals,sis),
172  send = refinementEnd(refinement,intervals,sis);
173  sit != send;
174  ++sit)
175  {
176  f.write(sit.coords(),*p);
177  // expand 2D-Vectors to 3D for VTK format
178  for(unsigned j = f.fieldInfo().size(); j < writecomps; j++)
179  p->write(0.0);
180  }
181  f.unbind();
182  }
183  }
184  }
185 
186 
187  protected:
189  virtual void countEntities(int &nvertices, int &ncells, int &ncorners);
190 
192  virtual void writeCellData(VTK::VTUWriter& writer);
193 
195  virtual void writeVertexData(VTK::VTUWriter& writer);
196 
198  virtual void writeGridPoints(VTK::VTUWriter& writer);
199 
201  virtual void writeGridCells(VTK::VTUWriter& writer);
202 
203  public:
204  using Base::addVertexData;
205  using Base::addCellData;
206 
207  private:
208  // hide addVertexData -- adding raw data directly without a VTKFunction
209  // currently does not make sense for subsampled meshes, as the higher order
210  // information is missing. See FS#676.
211  template<class V>
212  void addVertexData (const V& v, const std::string &name, int ncomps=1);
213  template<class V>
214  void addCellData (const V& v, const std::string &name, int ncomps=1);
215 
216  Dune::RefinementIntervals intervals;
217  bool coerceToSimplex;
218  };
219 
221  template <class GridView>
222  void SubsamplingVTKWriter<GridView>::countEntities(int &nvertices, int &ncells, int &ncorners)
223  {
224  nvertices = 0;
225  ncells = 0;
226  ncorners = 0;
227  for (CellIterator it=this->cellBegin(); it!=cellEnd(); ++it)
228  {
229  Refinement &refinement = buildRefinement<dim, ctype>(it->type(), subsampledGeometryType(it->type()));
230 
231  ncells += refinement.nElements(intervals);
232  nvertices += refinement.nVertices(intervals);
233  ncorners += refinement.nElements(intervals) * refinement.eBegin(intervals).vertexIndices().size();
234  }
235  }
236 
237 
239  template <class GridView>
241  {
242  if(celldata.size() == 0)
243  return;
244 
245  // Find the names of the first scalar and vector data fields.
246  // These will be marked as the default fields (the ones that ParaView shows
247  // when the file has just been opened).
248  std::string defaultScalarField, defaultVectorField;
249  std::tie(defaultScalarField, defaultVectorField) = this->getDataNames(celldata);
250 
251  writer.beginCellData(defaultScalarField, defaultVectorField);
252  writeData(writer,celldata,cellBegin(),cellEnd(),ncells,IteratorSelector<SubElementIterator>());
253  writer.endCellData();
254  }
255 
257  template <class GridView>
259  {
260  if(vertexdata.size() == 0)
261  return;
262 
263  // Find the names of the first scalar and vector data fields.
264  // These will be marked as the default fields (the ones that ParaView shows
265  // when the file has just been opened).
266  std::string defaultScalarField, defaultVectorField;
267  std::tie(defaultScalarField, defaultVectorField) = this->getDataNames(vertexdata);
268 
269  writer.beginPointData(defaultScalarField, defaultVectorField);
270  writeData(writer,vertexdata,cellBegin(),cellEnd(),nvertices,IteratorSelector<SubVertexIterator>());
271  writer.endPointData();
272  }
273 
275  template <class GridView>
277  {
278  writer.beginPoints();
279 
280  std::shared_ptr<VTK::DataArrayWriter> p
281  (writer.makeArrayWriter("Coordinates", 3, nvertices, this->coordPrecision()));
282  if(!p->writeIsNoop())
283  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
284  {
285  Refinement &refinement =
286  buildRefinement<dim, ctype>(i->type(),
287  subsampledGeometryType(i->type()));
288  for(SubVertexIterator sit = refinement.vBegin(intervals),
289  send = refinement.vEnd(intervals);
290  sit != send; ++sit)
291  {
292  FieldVector<ctype, dimw> coords = i->geometry().global(sit.coords());
293  for (int j=0; j<std::min(int(dimw),3); j++)
294  p->write(coords[j]);
295  for (int j=std::min(int(dimw),3); j<3; j++)
296  p->write(0.0);
297  }
298  }
299  // free the VTK::DataArrayWriter before touching the stream
300  p.reset();
301 
302  writer.endPoints();
303  }
304 
306  template <class GridView>
308  {
309  writer.beginCells();
310 
311  // connectivity
312  {
313  std::shared_ptr<VTK::DataArrayWriter> p1
314  (writer.makeArrayWriter("connectivity", 1, ncorners, VTK::Precision::int32));
315  // The offset within the index numbering
316  if(!p1->writeIsNoop()) {
317  int offset = 0;
318  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
319  {
320  GeometryType coercedToType = subsampledGeometryType(i->type());
321  Refinement &refinement =
322  buildRefinement<dim, ctype>(i->type(), coercedToType);
323  for(SubElementIterator sit = refinement.eBegin(intervals),
324  send = refinement.eEnd(intervals);
325  sit != send; ++sit)
326  {
327  IndexVector indices = sit.vertexIndices();
328  for(unsigned int ii = 0; ii < indices.size(); ++ii)
329  p1->write(offset+indices[VTK::renumber(coercedToType, ii)]);
330  }
331  offset += refinement.nVertices(intervals);
332  }
333  }
334  }
335 
336  // offsets
337  {
338  std::shared_ptr<VTK::DataArrayWriter> p2
339  (writer.makeArrayWriter("offsets", 1, ncells, VTK::Precision::int32));
340  if(!p2->writeIsNoop()) {
341  // The offset into the connectivity array
342  int offset = 0;
343  for (CellIterator i=cellBegin(); i!=cellEnd(); ++i)
344  {
345  Refinement &refinement =
346  buildRefinement<dim, ctype>(i->type(),
347  subsampledGeometryType(i->type()));
348  unsigned int verticesPerCell =
349  refinement.eBegin(intervals).vertexIndices().size();
350  for(int element = 0; element < refinement.nElements(intervals);
351  ++element)
352  {
353  offset += verticesPerCell;
354  p2->write(offset);
355  }
356  }
357  }
358  }
359 
360  // types
361  if (dim>1)
362  {
363  std::shared_ptr<VTK::DataArrayWriter> p3
364  (writer.makeArrayWriter("types", 1, ncells, VTK::Precision::uint8));
365  if(!p3->writeIsNoop())
366  for (CellIterator it=cellBegin(); it!=cellEnd(); ++it)
367  {
368  GeometryType coerceTo = subsampledGeometryType(it->type());
369  Refinement &refinement =
370  buildRefinement<dim, ctype>(it->type(), coerceTo);
371  int vtktype = VTK::geometryType(coerceTo);
372  for(int i = 0; i < refinement.nElements(intervals); ++i)
373  p3->write(vtktype);
374  }
375  }
376 
377  writer.endCells();
378  }
379 }
380 
381 #endif // DUNE_SUBSAMPLINGVTKWRITER_HH
Dune::VTK::VTUWriter::beginCellData
void beginCellData(const std::string &scalars="", const std::string &vectors="")
start CellData section
Definition: vtuwriter.hh:203
Dune::VTKWriter::celldata
std::list< VTKLocalFunction > celldata
Definition: vtkwriter.hh:1537
Dune::SubsamplingVTKWriter::refinementBegin
SubElementIterator refinementBegin(const Refinement &refinement, Dune::RefinementIntervals intervals, IteratorSelector< SubElementIterator >)
Definition: subsamplingvtkwriter.hh:116
Dune::SubsamplingVTKWriter::writeVertexData
virtual void writeVertexData(VTK::VTUWriter &writer)
write vertex data
Definition: subsamplingvtkwriter.hh:258
Dune::SubsamplingVTKWriter::DUNE_DEPRECATED_MSG
DUNE_DEPRECATED_MSG("SubsampligVTKWriter(GV,int,bool) is deprecated, use SubsamplingVTKWriter(GV,Dune::refinement{Intervals|Levels}(int),bool)") explicit SubsamplingVTKWriter(const GridView &gridView
Construct a SubsamplingVTKWriter working on a specific GridView.
Dune::VTK::VTUWriter::endPoints
void endPoints()
finish section for the point coordinates
Definition: vtuwriter.hh:247
Dune::SubsamplingVTKWriter::refinementBegin
SubVertexIterator refinementBegin(const Refinement &refinement, Dune::RefinementIntervals intervals, IteratorSelector< SubVertexIterator >)
Definition: subsamplingvtkwriter.hh:121
Dune::SubsamplingVTKWriter::writeData
void writeData(VTK::VTUWriter &writer, const Data &data, const Iterator begin, const Iterator end, int nentries, IteratorSelector< SubIterator > sis)
Definition: subsamplingvtkwriter.hh:137
Dune::VTK::VTUWriter
Dump a .vtu/.vtp files contents to a stream.
Definition: vtuwriter.hh:96
Dune::VTKWriter::cellBegin
CellIterator cellBegin() const
Definition: vtkwriter.hh:393
Dune::VTKWriter::nvertices
int nvertices
Definition: vtkwriter.hh:1545
Dune::VTK::VTUWriter::makeArrayWriter
DataArrayWriter * makeArrayWriter(const std::string &name, unsigned ncomps, unsigned nitems, Precision prec)
acquire a DataArrayWriter
Definition: vtuwriter.hh:378
Dune::VTK::FieldInfo::Type::vector
@ vector
vector-valued field (always 3D, will be padded if necessary)
Dune::VTKWriter::CellIterator
Iterator over the grids elements.
Definition: vtkwriter.hh:380
Dune::SubsamplingVTKWriter::countEntities
virtual void countEntities(int &nvertices, int &ncells, int &ncorners)
count the vertices, cells and corners
Definition: subsamplingvtkwriter.hh:222
Dune::VTKWriter::ncells
int ncells
Definition: vtkwriter.hh:1544
Dune::VTKWriter::cellEnd
CellIterator cellEnd() const
Definition: vtkwriter.hh:398
Dune::VTK::VTUWriter::beginPointData
void beginPointData(const std::string &scalars="", const std::string &vectors="")
start PointData section
Definition: vtuwriter.hh:165
Dune::VTKWriter::vertexdata
std::list< VTKLocalFunction > vertexdata
Definition: vtkwriter.hh:1538
Dune::SubsamplingVTKWriter
Writer for the output of subsampled grid functions in the vtk format.
Definition: subsamplingvtkwriter.hh:36
Dune::VTK::DataArrayWriter::writeIsNoop
virtual bool writeIsNoop() const
whether calls to write may be skipped
Definition: dataarraywriter.hh:87
Dune::VTK::DataArrayWriter::write
void write(T data)
write one element of data
Definition: dataarraywriter.hh:67
Dune::VTK::VTUWriter::beginPoints
void beginPoints()
start section for the point coordinates
Definition: vtuwriter.hh:236
vtkwriter.hh
Provides file i/o for the visualization toolkit.
Dune::VTK::VTUWriter::endCellData
void endCellData()
finish CellData section
Definition: vtuwriter.hh:218
Dune::SubsamplingVTKWriter::writeCellData
virtual void writeCellData(VTK::VTUWriter &writer)
write cell data
Definition: subsamplingvtkwriter.hh:240
Dune::VTK::nonconforming
@ nonconforming
Output non-conforming data.
Definition: common.hh:80
Dune::VTKWriter::ncorners
int ncorners
Definition: vtkwriter.hh:1546
Dune::GridView::dimension
@ dimension
The dimension of the grid.
Definition: common/gridview.hh:127
Dune::VTKWriter::FunctionIterator
std::list< VTKLocalFunction >::const_iterator FunctionIterator
Definition: vtkwriter.hh:372
Dune::VTK::FieldInfo::size
std::size_t size() const
The number of components in the data field.
Definition: common.hh:412
Dune::SubsamplingVTKWriter::writeGridPoints
virtual void writeGridPoints(VTK::VTUWriter &writer)
write the positions of vertices
Definition: subsamplingvtkwriter.hh:276
Dune::SubsamplingVTKWriter::refinementEnd
SubVertexIterator refinementEnd(const Refinement &refinement, Dune::RefinementIntervals intervals, IteratorSelector< SubVertexIterator >)
Definition: subsamplingvtkwriter.hh:131
Dune::VTK::FieldInfo::precision
Precision precision() const
The precision used for the output of the data field.
Definition: common.hh:418
Dune::VTKWriter::vertexEnd
VertexIterator vertexEnd() const
Definition: vtkwriter.hh:511
Dune::VTK::Precision::uint8
@ uint8
Dune::SubsamplingVTKWriter::SubsamplingVTKWriter
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
Dune::SubsamplingVTKWriter::writeGridCells
virtual void writeGridCells(VTK::VTUWriter &writer)
write the connectivity array
Definition: subsamplingvtkwriter.hh:307
Dune::VTK::VTUWriter::endCells
void endCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:283
Dune::VTK::geometryType
GeometryType geometryType(const Dune::GeometryType &t)
mapping from GeometryType to VTKGeometryType
Definition: common.hh:199
Dune::VTK::FieldInfo::type
Type type() const
The type of the data field.
Definition: common.hh:406
Dune::VTK::VTUWriter::endPointData
void endPointData()
finish PointData section
Definition: vtuwriter.hh:180
Dune::VTKWriter::outputtype
VTK::OutputType outputtype
Definition: vtkwriter.hh:1562
Dune::GridView::dimensionworld
@ dimensionworld
The dimension of the world the grid lives in.
Definition: common/gridview.hh:131
Dune::VTKWriter::vertexBegin
VertexIterator vertexBegin() const
Definition: vtkwriter.hh:504
Dune::GridView
Grid view abstract base class.
Definition: common/gridview.hh:59
Dune::SubsamplingVTKWriter::coerceToSimplex_
int bool coerceToSimplex_
Definition: subsamplingvtkwriter.hh:102
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:709
Dune::SubsamplingVTKWriter::level_
int level_
Definition: subsamplingvtkwriter.hh:102
Dune::VTK::FieldInfo::Type::tensor
@ tensor
tensor field (always 3x3)
Dune::SubsamplingVTKWriter::refinementEnd
SubElementIterator refinementEnd(const Refinement &refinement, Dune::RefinementIntervals intervals, IteratorSelector< SubElementIterator >)
Definition: subsamplingvtkwriter.hh:126
Dune::VTKWriter
Writer for the ouput of grid functions in the vtk format.
Definition: vtkwriter.hh:91
Dune::VTK::Precision
Precision
which precision to use when writing out data to vtk files
Definition: common.hh:319
Dune::VTK::Precision::float32
@ float32
Dune::VTK::renumber
int renumber(const Dune::GeometryType &t, int i)
renumber VTK <-> Dune
Definition: common.hh:234
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:645
Dune::VTK::FieldInfo
Descriptor struct for VTK fields.
Definition: common.hh:375
Dune::VTK::VTUWriter::beginCells
void beginCells()
start section for the grid cells/PolyData lines
Definition: vtuwriter.hh:272
vtuwriter.hh
Dune::VTK::Precision::int32
@ int32
Dune
Include standard header files.
Definition: agrid.hh:58
Dune::VTK::GeometryType
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:180
Dune::VTKWriter::coordPrecision
VTK::Precision coordPrecision() const
get the precision with which coordinates are written out
Definition: vtkwriter.hh:778
Dune::VTK::FieldInfo::Type::scalar
@ scalar