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