dune-fem 2.8.0
Loading...
Searching...
No Matches
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>
13
15
16namespace Dune
17{
18
19 namespace Fem
20 {
21
22 // Internal Forward Declarations
23 // -----------------------------
24
25 template< class GridPart, bool subsampling = false >
26 class VTKIO;
27
28
29
30 // VTKFunctionWrapper
31 // ------------------
32
33 template< class DF >
35 : public VTKFunction< typename DF::GridPartType::GridViewType >
36 {
38
39 public:
42
44 typedef typename DiscreteFunctionType::FunctionSpaceType FunctionSpaceType;
45
46 static const int dimRange = FunctionSpaceType::dimRange;
47 static const int dimDomain = FunctionSpaceType::dimDomain;
48
49 typedef typename FunctionSpaceType::DomainType DomainType;
50 typedef typename FunctionSpaceType::RangeType RangeType;
51
52 typedef typename DiscreteFunctionType::GridPartType GridPartType;
53 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
54
55 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
56
59 const std::string& dataName,
60 int component, bool vector, TypeOfField typeOfField )
61 : localFunction_( df ),
62 name_( ( dataName.size() > 0 ) ? dataName : df.name() ),
63 vector_( vector ),
64 typeOfField_( typeOfField ),
65 component_( component )
66 {}
67
70 {}
71
73 virtual int ncomps () const
74 {
75 return (!vector_) ? 1 : 3; // dimDomain;
76 }
77
80 virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
81 {
82 localFunction_.bind( e );
83 typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
84 RangeType val;
85 localFunction_.evaluate(xi,val);
86 localFunction_.unbind();
87
88 RangeFieldType outVal( 0 );
89 if (vector_)
90 {
91 if( comp <= dimDomain )
92 outVal = val[ comp + component_ ] ;
93 }
94 else
95 outVal = val[ component_ ] ;
96 if (typeOfField_ == TypeOfField::real || typeOfField_ == TypeOfField::complex_real )
97 return std::real( outVal );
98 else
99 return std::imag( outVal );
100 }
101
103 virtual std::string name () const
104 {
105 std::stringstream ret;
106 ret << name_;
107 if (typeOfField_ == TypeOfField::complex_real)
108 ret << "_real_";
109 if (typeOfField_ == TypeOfField::complex_imag)
110 ret << "_imag_";
111 if (vector_)
112 ret << "_vec" << component_;
113 else
114 ret << component_;
115 return ret.str();
116 }
117
118 private:
119 mutable LocalFunctionType localFunction_;
120 const std::string name_ ;
121 const bool vector_;
122 const TypeOfField typeOfField_;
123 const int component_;
124 };
125
126
127
128 // VTKIOBase
129 // ---------
130
132 template< class GridPart, bool subsampling >
134 {
136
137 protected:
138 typedef typename GridPart::GridViewType GridViewType;
139
140 class VTKWriter;
142
143 typedef typename std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
145
146 public:
147 typedef GridPart GridPartType;
148
149 typedef typename GridPartType::GridType GridType;
150 typedef typename GridPartType::IndexSetType IndexSetType;
151
152 protected:
154 : public VTKFunction< GridViewType >
155 {
157
158 public:
159 typedef typename GridViewType :: template Codim< 0 >::Entity EntityType;
160 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
161
163
166 const std::string& name,
167 const int rank, const int nThreads )
168 : iterators_( gridPart ), name_( name ), rank_( rank ), nThreads_( nThreads ) {}
169
171 virtual ~PartitioningData () {}
172
174 virtual int ncomps () const { return 1; }
175
178 virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
179 {
180 const int thread = iterators_.thread( e );
181 return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
182 }
183
185 virtual std::string name () const
186 {
187 return name_;
188 }
189
190 private:
191 ThreadIteratorType iterators_;
192 const std::string name_;
193 const int rank_;
194 const int nThreads_;
195
196 };
197
199 : public VTKFunction< GridViewType >
200 {
202
203 public:
204 typedef typename GridViewType :: template Codim< 0 >::Entity EntityType;
205 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
206
208
211
213 virtual ~VolumeData () {}
214
216 virtual int ncomps () const { return 1; }
217
220 virtual double evaluate ( int comp, const EntityType &e, const LocalCoordinateType &xi ) const
221 {
222 return e.geometry().volume();
223 }
224
226 virtual std::string name () const
227 {
228 return std::string("volume");
229 }
230 };
231
233 {
234 // 0 = none, 1 = MPI ranks only, 2 = ranks + threads, 3 = like 1 and also threads only
235 const std::string names[] = { "none", "rank", "rank+thread", "rank/thread" };
236 return parameter.getEnum( "fem.io.partitioning", names, 0 );
237 }
238
239 protected :
241 : gridPart_( gridPart ),
242 vtkWriter_( vtkWriter ),
244 {
245 static const std::string typeTable[] = { "ascii", "base64", "appended-raw", "appended-base64" };
246 static const VTK::OutputType typeValue[] = { VTK::ascii, VTK::base64, VTK::appendedraw, VTK::appendedbase64 };
247 type_ = typeValue[ parameter.getEnum( "fem.io.vtk.type", typeTable, 2 ) ];
248 }
249
250 void addPartitionData( const int myRank = -1 )
251 {
252 if( addPartition_ > 0 )
253 {
254 std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
255 vtkWriter_->addCellData( volumePtr );
256
257 const int rank = ( myRank < 0 ) ? gridPart_.comm().rank() : myRank ;
258 const int nThreads = ( addPartition_ > 1 ) ? MPIManager::maxThreads() : 1 ;
259 if( addPartition_ <= 2 )
260 {
261 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_, "rank", rank, nThreads) );
262 vtkWriter_->addCellData( dataRankPtr );
263 }
264 else
265 {
266 // rank only visualization
267 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_, "rank", rank, -1) );
268 vtkWriter_->addCellData( dataRankPtr );
269 // thread only visualization
270 std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(gridPart_, "thread", 0, nThreads) );
271 vtkWriter_->addCellData( dataThreadPtr );
272 }
273 addPartition_ = 0 ;
274 }
275 }
276
277 template < class DF >
278 static bool notComplex()
279 {
280 typedef typename DF::RangeFieldType RangeFieldType;
281 typedef typename Dune::FieldTraits< RangeFieldType >::real_type RealType;
282 return ! std::is_same< typename std::remove_cv<RangeFieldType>::type, std::complex<RealType> >::value;
283 }
284
285 public:
287 {
288 delete vtkWriter_;
289 }
290
292 const GridPartType &gridPart () const
293 {
294 return gridPart_;
295 }
296
297 template< class DF >
298 void addCellData( DF &df , const std::string& dataName = "" )
299 {
300 static const int dimRange = DF::FunctionSpaceType::dimRange;
301 for( int i = 0;i < dimRange; ++i )
302 {
303 if ( notComplex<DF>() )
304 {
305 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
307 vtkWriter_->addCellData( ptr );
308 }
309 else
310 {
311 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
313 vtkWriter_->addCellData( ptrR );
314 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
316 vtkWriter_->addCellData( ptrI );
317 }
318 }
319 }
320
321 template< class DF >
322 void addVectorCellData( DF &df,
323 const std::string& dataName = "" ,
324 int startPoint = 0 )
325 {
326 if ( notComplex<DF>() )
327 {
328 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
330 vtkWriter_->addCellData( ptr );
331 }
332 else
333 {
334 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
336 vtkWriter_->addCellData( ptrR );
337 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
339 vtkWriter_->addCellData( ptrI );
340 }
341 }
342
343 template< class DF >
344 void addVertexData( DF &df, const std::string& dataName = "" )
345 {
346 static const int dimRange = DF::FunctionSpaceType::dimRange;
347 std::string name = ( dataName.size() > 0 ) ? dataName : df.name() ;
348 for( int i = 0;i < dimRange; ++i )
349 {
350 if ( notComplex<DF>() )
351 {
352 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
354 vtkWriter_->addVertexData( ptr );
355 }
356 else
357 {
358 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
360 vtkWriter_->addVertexData( ptrR );
361 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
363 vtkWriter_->addVertexData( ptrI );
364 }
365 }
366 }
367
368 template< class DF >
370 const std::string& dataName = "" ,
371 int startPoint = 0 )
372 {
373 if ( notComplex<DF>() )
374 {
375 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
377 vtkWriter_->addVertexData( ptr );
378 }
379 else
380 {
381 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
383 vtkWriter_->addVertexData( ptrR );
384 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
386 vtkWriter_->addVertexData( ptrI );
387 }
388 }
389
390 void clear ()
391 {
392 vtkWriter_->clear();
393 }
394
395 std::string write ( const std::string &name, VTK::OutputType type )
396 {
398 size_t pos = name.find_last_of( '/' );
399 if( pos != name.npos )
400 return vtkWriter_->pwrite( name.substr( pos+1, name.npos ), name.substr( 0, pos ), "", type );
401 else
402 return vtkWriter_->write( name, type );
403 }
404
405 std::string write ( const std::string &name )
406 {
407 return write( name, type_ );
408 }
409
410 std::string pwrite ( const std::string &name,
411 const std::string &path,
412 const std::string &extendpath,
413 VTK::OutputType type )
414 {
416 return vtkWriter_->pwrite( name, path, extendpath, type );
417 }
418
419 std::string pwrite ( const std::string &name,
420 const std::string &path,
421 const std::string &extendpath )
422 {
423 return pwrite( name, path, extendpath, type_ );
424 }
425
426 std::string write ( const std::string &name,
427 VTK::OutputType type,
428 const int rank,
429 const int size )
430 {
431 addPartitionData( rank );
432 return vtkWriter_->write( name, type, rank, size );
433 }
434
435 std::string write ( const std::string &name,
436 const int rank,
437 const int size )
438 {
439 return write( name, type_, rank, size );
440 }
441
442 protected:
446 VTK::OutputType type_;
447 };
448
449
450
451 // VTKIOBase::VTKWriter
452 // --------------------
453
454 template< class GridPart, bool subsampling >
455 class VTKIOBase< GridPart, subsampling >::VTKWriter
456 : public Dune::VTKWriter< GridViewType >
457 {
458 typedef Dune::VTKWriter< GridViewType > BaseType;
459
460 public:
461 // make all write methods public for data convert
462 using BaseType::write;
463 using BaseType::pwrite;
464
467 VTK::DataMode dm = VTK::conforming )
468 : BaseType( gridPart.gridView(), dm )
469 {}
470 };
471
472
473
474 // VTKIOBase::SubSamplingVTKWriter
475 // -------------------------------
476
477 template< class GridPart, bool subsampling >
478 class VTKIOBase< GridPart, subsampling >::SubsamplingVTKWriter
479 : public Dune::SubsamplingVTKWriter< GridViewType >
480 {
481 typedef Dune::SubsamplingVTKWriter< GridViewType > BaseType;
482
483 public:
484 // make all write methods public for data convert
485 using BaseType::write;
486 using BaseType::pwrite;
487
490 Dune::RefinementIntervals intervals,
491 bool coerceToSimplex = false )
492 : BaseType( gridPart.gridView(), intervals, 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 {
507
508 typedef typename BaseType::VTKWriterType VTKWriterType;
509
510 public:
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 {
533
534 typedef typename BaseType::VTKWriterType VTKWriterType;
535
536 public:
538
539 explicit VTKIO ( const GridPartType &gridPart, Dune::RefinementIntervals intervals, bool coerceToSimplex, const ParameterReader &parameter = Parameter::container() )
540 : BaseType( gridPart, new VTKWriterType( gridPart, intervals, coerceToSimplex ), parameter )
541 {}
542
543 explicit VTKIO ( const GridPartType &gridPart, unsigned int level, bool coerceToSimplex, const ParameterReader &parameter = Parameter::container() )
544 : VTKIO( gridPart, Dune::refinementLevels(level), coerceToSimplex, parameter )
545 {}
546
547 VTKIO ( const GridPartType &gridPart, unsigned int level, const ParameterReader &parameter = Parameter::container() )
548 : VTKIO( gridPart, level, false, parameter )
549 {}
550
551 VTKIO ( const GridPartType &gridPart, int level, const ParameterReader &parameter = Parameter::container() ) DUNE_DEPRECATED_MSG( "pass level as unsigned int" )
552 : VTKIO( gridPart, level, false, parameter )
553 {}
554
555 VTKIO ( const GridPartType &gridPart, const ParameterReader &parameter = Parameter::container() )
556 : VTKIO( gridPart, 0, false, parameter )
557 {}
558
559 VTKIO ( const GridPartType &gridPart, bool coerceToSimplex, const ParameterReader &parameter = Parameter::container() )
560 : VTKIO( gridPart, 0, coerceToSimplex, parameter )
561 {}
562 };
563
564
565
566 // SubsamplingVTKIO
567 // ----------------
568
569 template< class GridPart >
571
572 } // namespace Fem
573
574} // namespace Dune
575
576#endif // #ifndef DUNE_FEM_VTKIO_HH
std::string path
Definition: readioparams.cc:156
double imag(const complex< Dune::Fem::Double > &x)
Definition: double.hh:995
double real(const complex< Dune::Fem::Double > &x)
Definition: double.hh:983
Definition: bindguard.hh:11
typename Impl::ConstLocalFunction< GridFunction >::Type ConstLocalFunction
Definition: const.hh:604
BasicParameterReader< std::function< const std::string *(const std::string &, const std::string *) > > ParameterReader
Definition: reader.hh:315
Definition: vtkio.hh:26
Definition: vtkio.hh:36
static const int dimDomain
Definition: vtkio.hh:47
TypeOfField
Definition: vtkio.hh:40
@ complex_real
Definition: vtkio.hh:40
@ real
Definition: vtkio.hh:40
@ complex_imag
Definition: vtkio.hh:40
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:55
static const int dimRange
Definition: vtkio.hh:46
virtual std::string name() const
get name
Definition: vtkio.hh:103
VTKFunctionWrapper(const DiscreteFunctionType &df, const std::string &dataName, int component, bool vector, TypeOfField typeOfField)
constructor taking discrete function
Definition: vtkio.hh:58
FunctionSpaceType::DomainType DomainType
Definition: vtkio.hh:49
virtual ~VTKFunctionWrapper()
virtual destructor
Definition: vtkio.hh:69
DiscreteFunctionType::GridPartType GridPartType
Definition: vtkio.hh:52
virtual int ncomps() const
return number of components
Definition: vtkio.hh:73
FunctionSpaceType::RangeType RangeType
Definition: vtkio.hh:50
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: vtkio.hh:53
DF DiscreteFunctionType
Definition: vtkio.hh:41
ConstLocalFunction< DF > LocalFunctionType
Definition: vtkio.hh:43
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:80
DiscreteFunctionType::FunctionSpaceType FunctionSpaceType
Definition: vtkio.hh:44
Output using VTK.
Definition: vtkio.hh:134
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath, VTK::OutputType type)
Definition: vtkio.hh:410
int addPartition_
Definition: vtkio.hh:445
VTKIOBase(const GridPartType &gridPart, VTKWriterType *vtkWriter, const ParameterReader &parameter=Parameter::container())
Definition: vtkio.hh:240
std::string write(const std::string &name, VTK::OutputType type)
Definition: vtkio.hh:395
std::string write(const std::string &name, VTK::OutputType type, const int rank, const int size)
Definition: vtkio.hh:426
void addPartitionData(const int myRank=-1)
Definition: vtkio.hh:250
std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type VTKWriterType
Definition: vtkio.hh:144
VTK::OutputType type_
Definition: vtkio.hh:446
std::string pwrite(const std::string &name, const std::string &path, const std::string &extendpath)
Definition: vtkio.hh:419
static bool notComplex()
Definition: vtkio.hh:278
GridPart GridPartType
Definition: vtkio.hh:147
GridPartType::GridType GridType
Definition: vtkio.hh:149
void addVertexData(DF &df, const std::string &dataName="")
Definition: vtkio.hh:344
void addVectorCellData(DF &df, const std::string &dataName="", int startPoint=0)
Definition: vtkio.hh:322
const GridPartType & gridPart() const
return grid part
Definition: vtkio.hh:292
std::string write(const std::string &name)
Definition: vtkio.hh:405
std::string write(const std::string &name, const int rank, const int size)
Definition: vtkio.hh:435
GridPartType::IndexSetType IndexSetType
Definition: vtkio.hh:150
int getPartitionParameter(const ParameterReader &parameter=Parameter::container()) const
Definition: vtkio.hh:232
void clear()
Definition: vtkio.hh:390
~VTKIOBase()
Definition: vtkio.hh:286
const GridPartType & gridPart_
Definition: vtkio.hh:443
GridPart::GridViewType GridViewType
Definition: vtkio.hh:138
VTKWriterType * vtkWriter_
Definition: vtkio.hh:444
void addVectorVertexData(DF &df, const std::string &dataName="", int startPoint=0)
Definition: vtkio.hh:369
void addCellData(DF &df, const std::string &dataName="")
Definition: vtkio.hh:298
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:178
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:159
virtual ~PartitioningData()
virtual destructor
Definition: vtkio.hh:171
virtual int ncomps() const
return number of components
Definition: vtkio.hh:174
virtual std::string name() const
get name
Definition: vtkio.hh:185
PartitioningData(const GridPartType &gridPart, const std::string &name, const int rank, const int nThreads)
constructor taking discrete function
Definition: vtkio.hh:165
DomainDecomposedIteratorStorage< GridPartType > ThreadIteratorType
Definition: vtkio.hh:162
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:160
Definition: vtkio.hh:200
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:204
VolumeData()
constructor taking discrete function
Definition: vtkio.hh:210
virtual std::string name() const
get name
Definition: vtkio.hh:226
virtual int ncomps() const
return number of components
Definition: vtkio.hh:216
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:220
DomainDecomposedIteratorStorage< GridPartType > ThreadIteratorType
Definition: vtkio.hh:207
virtual ~VolumeData()
virtual destructor
Definition: vtkio.hh:213
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Definition: vtkio.hh:205
Definition: vtkio.hh:457
VTKWriter(const GridPartType &gridPart, VTK::DataMode dm=VTK::conforming)
constructor
Definition: vtkio.hh:466
SubsamplingVTKWriter(const GridPartType &gridPart, Dune::RefinementIntervals intervals, bool coerceToSimplex=false)
constructor
Definition: vtkio.hh:489
BaseType::GridPartType GridPartType
Definition: vtkio.hh:511
Definition: vtkio.hh:530
VTKIO(const GridPartType &gridPart, Dune::RefinementIntervals intervals, bool coerceToSimplex, const ParameterReader &parameter=Parameter::container())
Definition: vtkio.hh:539
BaseType::GridPartType GridPartType
Definition: vtkio.hh:537
static ParameterContainer & container()
Definition: io/parameter.hh:193
static int maxThreads()
return maximal number of threads possible in the current run
Definition: mpimanager.hh:418
int thread(const EntityType &entity) const
return thread number this entity belongs to
Definition: threaditeratorstorage.hh:125