dune-fem  2.4.1-rc
vtkio.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_VTKIO_HH
2 #define DUNE_FEM_VTKIO_HH
3 
4 #include <dune/common/typetraits.hh>
5 
6 #include <dune/common/deprecated.hh>
7 
8 #include <dune/grid/io/file/vtk/vtkwriter.hh>
9 #include <dune/grid/io/file/vtk/subsamplingvtkwriter.hh>
10 
11 #include <dune/fem/version.hh>
12 #include <dune/fem/io/parameter.hh>
13 
16 
17 namespace Dune
18 {
19 
20  namespace Fem
21  {
22 
23  // Internal Forward Declarations
24  // -----------------------------
25 
26  template< class GridPart, bool subsampling = false >
27  class VTKIO;
28 
29 
30 
31  // VTKFunctionWrapper
32  // ------------------
33 
34  template< class DF >
36  : public VTKFunction< typename DF::GridPartType::GridViewType >
37  {
39 
40  public:
43 
44  typedef typename DiscreteFunctionType::LocalFunctionType LocalFunctionType;
45  typedef typename DiscreteFunctionType::FunctionSpaceType FunctionSpaceType;
46 
47  static const int dimRange = FunctionSpaceType::dimRange;
48  static const int dimDomain = FunctionSpaceType::dimDomain;
49 
50  typedef typename FunctionSpaceType::DomainType DomainType;
51  typedef typename FunctionSpaceType::RangeType RangeType;
52 
53  typedef typename DiscreteFunctionType::GridPartType GridPartType;
54  typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
55 
56  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
57 
59  VTKFunctionWrapper ( const DiscreteFunctionType& df,
60  const std::string& dataName,
61  int component, bool vector, TypeOfField typeOfField )
62  : discFunc_( df ),
63  name_( ( dataName.size() > 0 ) ? dataName : df.name() ),
64  vector_( vector ),
65  typeOfField_( typeOfField ),
66  component_( component )
67  {}
68 
71  {}
72 
74  virtual int ncomps () const
75  {
76  return (!vector_) ? 1 : 3; // dimDomain;
77  }
78 
81  virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
82  {
83  typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
84  const LocalFunctionType lf = discFunc_.localFunction(e);
85  RangeType val;
86  lf.evaluate(xi,val);
87  RangeFieldType outVal( 0 );
88  if (vector_)
89  {
90  if( comp <= dimDomain )
91  outVal = val[ comp + component_ ] ;
92  }
93  else
94  outVal = val[ component_ ] ;
95  if (typeOfField_ == TypeOfField::real || typeOfField_ == TypeOfField::complex_real )
96  return std::real( outVal );
97  else
98  return std::imag( outVal );
99  }
100 
102  virtual std::string name () const
103  {
104  std::stringstream ret;
105  ret << name_;
106  if (typeOfField_ == TypeOfField::complex_real)
107  ret << "_real_";
108  if (typeOfField_ == TypeOfField::complex_imag)
109  ret << "_imag_";
110  if (vector_)
111  ret << "_vec" << component_;
112  else
113  ret << component_;
114  return ret.str();
115  }
116 
117  private:
118  const DiscreteFunctionType &discFunc_;
119  const std::string name_ ;
120  const bool vector_;
121  const TypeOfField typeOfField_;
122  const int component_;
123  };
124 
125 
126 
127  // VTKIOBase
128  // ---------
129 
131  template< class GridPart, bool subsampling >
132  class VTKIOBase
133  {
135 
136  protected:
137  typedef typename GridPart::GridViewType GridViewType;
138 
139  class VTKWriter;
140  class SubsamplingVTKWriter;
141 
142  typedef typename conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
144 
145  public:
146  typedef GridPart GridPartType;
147 
148  typedef typename GridPartType::GridType GridType;
149  typedef typename GridPartType::IndexSetType IndexSetType;
150 
151  protected:
153  : public VTKFunction< GridViewType >
154  {
155  typedef PartitioningData ThisType;
156 
157  public:
158  typedef typename GridViewType :: template Codim< 0 >::Entity EntityType;
159  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
160 
162 
164  PartitioningData( const GridPartType& gridPart,
165  const std::string& name,
166  const int rank, const int nThreads )
167  : iterators_( gridPart ), name_( name ), rank_( rank ), nThreads_( nThreads ) {}
168 
170  virtual ~PartitioningData () {}
171 
173  virtual int ncomps () const { return 1; }
174 
177  virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
178  {
179  const int thread = iterators_.thread( e );
180  return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
181  }
182 
184  virtual std::string name () const
185  {
186  return name_;
187  }
188 
189  private:
190  ThreadIteratorType iterators_;
191  const std::string name_;
192  const int rank_;
193  const int nThreads_;
194 
195  };
196 
198  : public VTKFunction< GridViewType >
199  {
200  typedef PartitioningData ThisType;
201 
202  public:
203  typedef typename GridViewType :: template Codim< 0 >::Entity EntityType;
204  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
205 
207 
210 
212  virtual ~VolumeData () {}
213 
215  virtual int ncomps () const { return 1; }
216 
219  virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
220  {
221  return e.geometry().volume();
222  }
223 
225  virtual std::string name () const
226  {
227  return std::string("volume");
228  }
229  };
230 
231  int getPartitionParameter ( const ParameterReader &parameter = Parameter::container() ) const
232  {
233  // 0 = none, 1 = MPI ranks only, 2 = ranks + threads, 3 = like 1 and also threads only
234  const std::string names[] = { "none", "rank", "rank+thread", "rank/thread" };
235  return parameter.getEnum( "fem.io.partitioning", names, 0 );
236  }
237 
238  protected :
239  VTKIOBase ( const GridPartType &gridPart, VTKWriterType *vtkWriter, const ParameterReader &parameter = Parameter::container() )
240  : gridPart_( gridPart ),
241  vtkWriter_( vtkWriter ),
242  addPartition_( getPartitionParameter( parameter ) )
243  {
244  static const std::string typeTable[] = { "ascii", "base64", "appended-raw", "appended-base64" };
245  static const VTK::OutputType typeValue[] = { VTK::ascii, VTK::base64, VTK::appendedraw, VTK::appendedbase64 };
246  type_ = typeValue[ parameter.getEnum( "fem.io.vtk.type", typeTable, 2 ) ];
247  }
248 
249  void addPartitionData( const int myRank = -1 )
250  {
251  if( addPartition_ > 0 )
252  {
253  std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
254  vtkWriter_->addCellData( volumePtr );
255 
256  const int rank = ( myRank < 0 ) ? gridPart_.comm().rank() : myRank ;
257  const int nThreads = ( addPartition_ > 1 ) ? ThreadManager::maxThreads() : 1 ;
258  if( addPartition_ <= 2 )
259  {
260  std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_, "rank", rank, nThreads) );
261  vtkWriter_->addCellData( dataRankPtr );
262  }
263  else
264  {
265  // rank only visualization
266  std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_, "rank", rank, -1) );
267  vtkWriter_->addCellData( dataRankPtr );
268  // thread only visualization
269  std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(gridPart_, "thread", 0, nThreads) );
270  vtkWriter_->addCellData( dataThreadPtr );
271  }
272  addPartition_ = 0 ;
273  }
274  }
275 
276  template < class DF >
277  static bool notComplex()
278  {
279  typedef typename DF::RangeFieldType RangeFieldType;
280  typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
281  return ! std::is_same< typename std::remove_cv<RangeFieldType>::type, std::complex<RealType> >::value;
282  }
283 
284  public:
286  {
287  delete vtkWriter_;
288  }
289 
291  const GridPartType &gridPart () const
292  {
293  return gridPart_;
294  }
295 
296  template< class DF >
297  void addCellData( DF &df , const std::string& dataName = "" )
298  {
299  static const int dimRange = DF::FunctionSpaceType::dimRange;
300  for( int i = 0;i < dimRange; ++i )
301  {
302  if ( notComplex<DF>() )
303  {
304  std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
306  vtkWriter_->addCellData( ptr );
307  }
308  else
309  {
310  std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
312  vtkWriter_->addCellData( ptrR );
313  std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
315  vtkWriter_->addCellData( ptrI );
316  }
317  }
318  }
319 
320  template< class DF >
321  void addVectorCellData( DF &df,
322  const std::string& dataName = "" ,
323  int startPoint = 0 )
324  {
325  if ( notComplex<DF>() )
326  {
327  std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
329  vtkWriter_->addCellData( ptr );
330  }
331  else
332  {
333  std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
335  vtkWriter_->addCellData( ptrR );
336  std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
338  vtkWriter_->addCellData( ptrI );
339  }
340  }
341 
342  template< class DF >
343  void addVertexData( DF &df, const std::string& dataName = "" )
344  {
345  static const int dimRange = DF::FunctionSpaceType::dimRange;
346  std::string name = ( dataName.size() > 0 ) ? dataName : df.name() ;
347  for( int i = 0;i < dimRange; ++i )
348  {
349  if ( notComplex<DF>() )
350  {
351  std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
353  vtkWriter_->addVertexData( ptr );
354  }
355  else
356  {
357  std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
359  vtkWriter_->addVertexData( ptrR );
360  std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
362  vtkWriter_->addVertexData( ptrI );
363  }
364  }
365  }
366 
367  template< class DF >
368  void addVectorVertexData( DF &df,
369  const std::string& dataName = "" ,
370  int startPoint = 0 )
371  {
372  if ( notComplex<DF>() )
373  {
374  std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
376  vtkWriter_->addVertexData( ptr );
377  }
378  else
379  {
380  std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
382  vtkWriter_->addVertexData( ptrR );
383  std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
385  vtkWriter_->addVertexData( ptrI );
386  }
387  }
388 
389  void clear ()
390  {
391  vtkWriter_->clear();
392  }
393 
394  std::string write ( const std::string &name, VTK::OutputType type )
395  {
396  addPartitionData();
397  size_t pos = name.find_last_of( '/' );
398  if( pos != name.npos )
399  return vtkWriter_->pwrite( name.substr( pos+1, name.npos ), name.substr( 0, pos ), "", type );
400  else
401  return vtkWriter_->write( name, type );
402  }
403 
404  std::string write ( const std::string &name )
405  {
406  return write( name, type_ );
407  }
408 
409  std::string pwrite ( const std::string &name,
410  const std::string &path,
411  const std::string &extendpath,
412  VTK::OutputType type )
413  {
414  addPartitionData();
415  return vtkWriter_->pwrite( name, path, extendpath, type );
416  }
417 
418  std::string pwrite ( const std::string &name,
419  const std::string &path,
420  const std::string &extendpath )
421  {
422  return pwrite( name, path, extendpath, type_ );
423  }
424 
425  std::string write ( const std::string &name,
426  VTK::OutputType type,
427  const int rank,
428  const int size )
429  {
430  addPartitionData( rank );
431  return vtkWriter_->write( name, type, rank, size );
432  }
433 
434  std::string write ( const std::string &name,
435  const int rank,
436  const int size )
437  {
438  return write( name, type_, rank, size );
439  }
440 
441  protected:
442  const GridPartType &gridPart_;
443  VTKWriterType *vtkWriter_;
445  VTK::OutputType type_;
446  };
447 
448 
449 
450  // VTKIOBase::VTKWriter
451  // --------------------
452 
453  template< class GridPart, bool subsampling >
454  class VTKIOBase< GridPart, subsampling >::VTKWriter
455  : public Dune::VTKWriter< GridViewType >
456  {
457  typedef Dune::VTKWriter< GridViewType > BaseType;
458 
459  public:
460  // make all write methods public for data convert
461  using BaseType::write;
462  using BaseType::pwrite;
463 
465  VTKWriter( const GridPartType &gridPart,
466  VTK::DataMode dm = VTK::conforming )
467  : BaseType( static_cast< GridViewType > ( gridPart ), dm )
468  {}
469  };
470 
471 
472 
473  // VTKIOBase::SubSamplingVTKWriter
474  // -------------------------------
475 
476  template< class GridPart, bool subsampling >
477  class VTKIOBase< GridPart, subsampling >::SubsamplingVTKWriter
478  : public Dune::SubsamplingVTKWriter< GridViewType >
479  {
480  typedef Dune::SubsamplingVTKWriter< GridViewType > BaseType;
481 
482  public:
483  // make all write methods public for data convert
484  using BaseType::write;
485  using BaseType::pwrite;
486 
488  SubsamplingVTKWriter( const GridPartType &gridPart,
489  const int level,
490  bool coerceToSimplex = false )
491  : BaseType( static_cast< GridViewType > ( gridPart ),
492  level, coerceToSimplex )
493  {}
494  };
495 
496 
497 
498  // VTKIO (without subsampling)
499  // ---------------------------
500 
501  template< class GridPart >
502  class VTKIO< GridPart, false >
503  : public VTKIOBase< GridPart, false >
504  {
505  typedef VTKIO< GridPart > ThisType;
507 
508  typedef typename BaseType::VTKWriterType VTKWriterType;
509 
510  public:
511  typedef typename BaseType::GridPartType GridPartType;
512 
513  VTKIO ( const GridPartType &gridPart, VTK::DataMode dm, const ParameterReader &parameter = Parameter::container() )
514  : BaseType( gridPart, new VTKWriterType( gridPart, dm ), parameter )
515  {}
516 
517  explicit VTKIO ( const GridPartType &gridPart, const ParameterReader &parameter = Parameter::container() )
518  : VTKIO( gridPart, VTK::conforming, parameter )
519  {}
520  };
521 
522 
523 
524  // VTKIO (with subsampling)
525  // ------------------------
526 
527  template< class GridPart >
528  class VTKIO< GridPart, true >
529  : public VTKIOBase< GridPart , true >
530  {
531  typedef VTKIO< GridPart, true > ThisType;
532  typedef VTKIOBase< GridPart, true > BaseType;
533 
534  typedef typename BaseType::VTKWriterType VTKWriterType;
535 
536  public:
537  typedef typename BaseType::GridPartType GridPartType;
538 
539  explicit VTKIO ( const GridPartType &gridPart, unsigned int level, bool coerceToSimplex, const ParameterReader &parameter = Parameter::container() )
540  : BaseType( gridPart, new VTKWriterType( gridPart, level, coerceToSimplex ), parameter )
541  {}
542 
543  VTKIO ( const GridPartType &gridPart, unsigned int level, const ParameterReader &parameter = Parameter::container() )
544  : VTKIO( gridPart, level, false, parameter )
545  {}
546 
547  VTKIO ( const GridPartType &gridPart, int level, const ParameterReader &parameter = Parameter::container() ) DUNE_DEPRECATED_MSG( "pass level as unsigned int" )
548  : VTKIO( gridPart, level, false, parameter )
549  {}
550 
551  VTKIO ( const GridPartType &gridPart, const ParameterReader &parameter = Parameter::container() )
552  : VTKIO( gridPart, 0, false, parameter )
553  {}
554 
555  VTKIO ( const GridPartType &gridPart, bool coerceToSimplex, const ParameterReader &parameter = Parameter::container() )
556  : VTKIO( gridPart, 0, coerceToSimplex, parameter )
557  {}
558 
559  };
560 
561 
562 
563  // SubsamplingVTKIO
564  // ----------------
565 
566  template< class GridPart >
568 
569  } // namespace Fem
570 
571 } // namespace Dune
572 
573 #endif // #ifndef DUNE_FEM_VTKIO_HH
GridPartType::IndexSetType IndexSetType
Definition: vtkio.hh:149
VTKWriterType * vtkWriter_
Definition: vtkio.hh:443
double imag(const std::complex< Double > &x)
Definition: double.hh:900
VTKIO(const GridPartType &gridPart, unsigned int level, bool coerceToSimplex, const ParameterReader &parameter=Parameter::container())
Definition: vtkio.hh:539
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:177
DiscreteFunctionType::GridPartType GridPartType
Definition: vtkio.hh:53
VTKIOBase(const GridPartType &gridPart, VTKWriterType *vtkWriter, const ParameterReader &parameter=Parameter::container())
Definition: vtkio.hh:239
virtual int ncomps() const
return number of components
Definition: vtkio.hh:173
GridPart GridPartType
Definition: vtkio.hh:146
void addPartitionData(const int myRank=-1)
Definition: vtkio.hh:249
DF DiscreteFunctionType
Definition: vtkio.hh:42
GridPartType::GridType GridType
Definition: vtkio.hh:148
DomainDecomposedIteratorStorage< GridPartType > ThreadIteratorType
Definition: vtkio.hh:161
Definition: vtkio.hh:35
void addVectorCellData(DF &df, const std::string &dataName="", int startPoint=0)
Definition: vtkio.hh:321
static const int dimDomain
Definition: vtkio.hh:48
BaseType::GridPartType GridPartType
Definition: vtkio.hh:511
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: vtkio.hh:54
Definition: vtkio.hh:528
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:219
Definition: vtkio.hh:454
virtual ~PartitioningData()
virtual destructor
Definition: vtkio.hh:170
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath, VTK::OutputType type)
Definition: vtkio.hh:409
conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type VTKWriterType
Definition: vtkio.hh:140
GridPart::GridViewType GridViewType
Definition: vtkio.hh:137
FunctionSpaceType::DomainType DomainType
Definition: vtkio.hh:50
virtual std::string name() const
get name
Definition: vtkio.hh:102
VTKFunctionWrapper(const DiscreteFunctionType &df, const std::string &dataName, int component, bool vector, TypeOfField typeOfField)
constructor taking discrete function
Definition: vtkio.hh:59
TypeOfField
Definition: vtkio.hh:41
VTKWriter(const GridPartType &gridPart, VTK::DataMode dm=VTK::conforming)
constructor
Definition: vtkio.hh:465
BaseType::GridPartType GridPartType
Definition: vtkio.hh:537
std::string write(const std::string &name)
Definition: vtkio.hh:404
void addVertexData(DF &df, const std::string &dataName="")
Definition: vtkio.hh:343
FunctionSpaceType::RangeType RangeType
Definition: vtkio.hh:51
VolumeData()
constructor taking discrete function
Definition: vtkio.hh:209
const GridPartType & gridPart_
Definition: vtkio.hh:442
Definition: coordinate.hh:4
DiscreteFunctionType::FunctionSpaceType FunctionSpaceType
Definition: vtkio.hh:45
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:159
virtual ~VolumeData()
virtual destructor
Definition: vtkio.hh:212
static ParameterContainer & container()
Definition: io/parameter.hh:190
std::string write(const std::string &name, VTK::OutputType type)
Definition: vtkio.hh:394
std::string write(const std::string &name, VTK::OutputType type, const int rank, const int size)
Definition: vtkio.hh:425
~VTKIOBase()
Definition: vtkio.hh:285
static bool notComplex()
Definition: vtkio.hh:277
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:158
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:81
int getPartitionParameter(const ParameterReader &parameter=Parameter::container()) const
Definition: vtkio.hh:231
std::string path
Definition: readioparams.cc:155
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:56
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:203
PartitioningData(const GridPartType &gridPart, const std::string &name, const int rank, const int nThreads)
constructor taking discrete function
Definition: vtkio.hh:164
static int maxThreads()
return maximal number of threads possbile in the current run
Definition: threadmanager.hh:202
void addCellData(DF &df, const std::string &dataName="")
Definition: vtkio.hh:297
double real(const std::complex< Double > &x)
Definition: double.hh:890
std::string write(const std::string &name, const int rank, const int size)
Definition: vtkio.hh:434
virtual int ncomps() const
return number of components
Definition: vtkio.hh:215
virtual int ncomps() const
return number of components
Definition: vtkio.hh:74
virtual std::string name() const
get name
Definition: vtkio.hh:184
Definition: vtkio.hh:197
const GridPartType & gridPart() const
return grid part
Definition: vtkio.hh:291
DiscreteFunctionType::LocalFunctionType LocalFunctionType
Definition: vtkio.hh:44
Definition: vtkio.hh:27
/brief Output using VTK
Definition: vtkio.hh:132
virtual std::string name() const
get name
Definition: vtkio.hh:225
static const int dimRange
Definition: vtkio.hh:47
int addPartition_
Definition: vtkio.hh:444
void addVectorVertexData(DF &df, const std::string &dataName="", int startPoint=0)
Definition: vtkio.hh:368
SubsamplingVTKWriter(const GridPartType &gridPart, const int level, bool coerceToSimplex=false)
constructor
Definition: vtkio.hh:488
VTK::OutputType type_
Definition: vtkio.hh:445
void clear()
Definition: vtkio.hh:389
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath)
Definition: vtkio.hh:418
DomainDecomposedIteratorStorage< GridPartType > ThreadIteratorType
Definition: vtkio.hh:206
virtual ~VTKFunctionWrapper()
virtual destructor
Definition: vtkio.hh:70
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:204