1#ifndef DUNE_FEM_VTKIO_HH
2#define DUNE_FEM_VTKIO_HH
11#include <dune/fem/version.hh>
12#include <dune/fem/io/parameter.hh>
14#include <dune/fem/misc/threads/domainthreaditerator.hh>
25 template<
class Gr
idPart,
bool subsampling = false >
34 class VTKFunctionWrapper
35 :
public VTKFunction< typename DF::GridPartType::GridViewType >
37 typedef VTKFunctionWrapper< DF > ThisType;
40 enum TypeOfField { real, complex_real, complex_imag };
43 typedef ConstLocalFunction< DF > LocalFunctionType;
53 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
55 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
59 const std::string& dataName,
60 int component,
bool vector, TypeOfField typeOfField )
61 : localFunction_( df ),
62 name_( ( dataName.
size() > 0 ) ? dataName : df.name() ),
64 typeOfField_( typeOfField ),
65 component_( component )
69 virtual ~VTKFunctionWrapper ()
73 virtual int ncomps ()
const
75 return (!vector_) ? 1 : 3;
80 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const
82 localFunction_.bind( e );
83 typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
85 localFunction_.evaluate(xi,val);
86 localFunction_.unbind();
88 RangeFieldType outVal( 0 );
91 if( comp <= dimDomain )
92 outVal = val[ comp + component_ ] ;
95 outVal = val[ component_ ] ;
96 if (typeOfField_ == TypeOfField::real || typeOfField_ == TypeOfField::complex_real )
97 return std::real( outVal );
99 return std::imag( outVal );
103 virtual std::string name ()
const
105 std::stringstream ret;
107 if (typeOfField_ == TypeOfField::complex_real)
109 if (typeOfField_ == TypeOfField::complex_imag)
112 ret <<
"_vec" << component_;
119 mutable LocalFunctionType localFunction_;
120 const std::string name_ ;
122 const TypeOfField typeOfField_;
123 const int component_;
132 template<
class Gr
idPart,
bool subsampling >
138 typedef typename GridPart::GridViewType GridViewType;
141 class SubsamplingVTKWriter;
143 typedef typename std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
147 typedef GridPart GridPartType;
149 typedef typename GridPartType::GridType GridType;
153 class PartitioningData
160 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
165 PartitioningData(
const GridPartType&
gridPart,
166 const std::string& name,
167 const int rank,
const int nThreads )
168 : iterators_(
gridPart ), name_( name ), rank_( rank ), nThreads_( nThreads ) {}
171 virtual ~PartitioningData () {}
174 virtual int ncomps ()
const {
return 1; }
178 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const
180 const int thread = iterators_.thread( e );
181 return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
185 virtual std::string name ()
const
191 ThreadIteratorType iterators_;
192 const std::string name_;
205 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
213 virtual ~VolumeData () {}
216 virtual int ncomps ()
const {
return 1; }
220 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const
222 return e.geometry().volume();
226 virtual std::string name ()
const
228 return std::string(
"volume");
232 int getPartitionParameter (
const ParameterReader ¶meter = Parameter::container() )
const
235 const std::string names[] = {
"none",
"rank",
"rank+thread",
"rank/thread" };
236 return parameter.getEnum(
"fem.io.partitioning", names, 0 );
240 VTKIOBase (
const GridPartType &
gridPart, VTKWriterType *vtkWriter,
const ParameterReader ¶meter = Parameter::container() )
242 vtkWriter_( vtkWriter ),
243 addPartition_( getPartitionParameter( parameter ) )
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 ) ];
250 void addPartitionData(
const int myRank = -1 )
252 if( addPartition_ > 0 )
254 std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
255 vtkWriter_->addCellData( volumePtr );
257 const int rank = ( myRank < 0 ) ? gridPart_.comm().rank() : myRank ;
258 const int nThreads = ( addPartition_ > 1 ) ? MPIManager::maxThreads() : 1 ;
259 if( addPartition_ <= 2 )
261 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_,
"rank", rank, nThreads) );
262 vtkWriter_->addCellData( dataRankPtr );
267 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_,
"rank", rank, -1) );
268 vtkWriter_->addCellData( dataRankPtr );
270 std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(gridPart_,
"thread", 0, nThreads) );
271 vtkWriter_->addCellData( dataThreadPtr );
277 template <
class DF >
278 static bool notComplex()
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;
298 void addCellData( DF &df ,
const std::string& dataName =
"" )
300 static const int dimRange = DF::FunctionSpaceType::dimRange;
301 for(
int i = 0;i < dimRange; ++i )
303 if ( notComplex<DF>() )
305 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
306 false, VTKFunctionWrapper<DF>::TypeOfField::real) );
307 vtkWriter_->addCellData( ptr );
311 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
312 false, VTKFunctionWrapper<DF>::TypeOfField::complex_real) );
313 vtkWriter_->addCellData( ptrR );
314 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
315 false, VTKFunctionWrapper<DF>::TypeOfField::complex_imag) );
316 vtkWriter_->addCellData( ptrI );
322 void addVectorCellData( DF &df,
323 const std::string& dataName =
"" ,
326 if ( notComplex<DF>() )
328 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
329 true, VTKFunctionWrapper<DF>::TypeOfField::real) );
330 vtkWriter_->addCellData( ptr );
334 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
335 true, VTKFunctionWrapper<DF>::TypeOfField::complex_real) );
336 vtkWriter_->addCellData( ptrR );
337 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
338 true, VTKFunctionWrapper<DF>::TypeOfField::complex_imag) );
339 vtkWriter_->addCellData( ptrI );
344 void addVertexData( DF &df,
const std::string& dataName =
"" )
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 )
350 if ( notComplex<DF>() )
352 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
353 false, VTKFunctionWrapper<DF>::TypeOfField::real) );
354 vtkWriter_->addVertexData( ptr );
358 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
359 false, VTKFunctionWrapper<DF>::TypeOfField::complex_real) );
360 vtkWriter_->addVertexData( ptrR );
361 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, i,
362 false, VTKFunctionWrapper<DF>::TypeOfField::complex_imag) );
363 vtkWriter_->addVertexData( ptrI );
369 void addVectorVertexData( DF &df,
370 const std::string& dataName =
"" ,
373 if ( notComplex<DF>() )
375 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
376 true, VTKFunctionWrapper<DF>::TypeOfField::real) );
377 vtkWriter_->addVertexData( ptr );
381 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
382 true, VTKFunctionWrapper<DF>::TypeOfField::complex_real) );
383 vtkWriter_->addVertexData( ptrR );
384 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<VTKFunctionWrapper< DF > >( df, dataName, startPoint,
385 true, VTKFunctionWrapper<DF>::TypeOfField::complex_imag) );
386 vtkWriter_->addVertexData( ptrI );
395 std::string write (
const std::string &name, VTK::OutputType type )
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 );
402 return vtkWriter_->write( name, type );
405 std::string write (
const std::string &name )
407 return write( name, type_ );
410 std::string pwrite (
const std::string &name,
411 const std::string &path,
412 const std::string &extendpath,
413 VTK::OutputType type )
416 return vtkWriter_->pwrite( name, path, extendpath, type );
419 std::string pwrite (
const std::string &name,
420 const std::string &path,
421 const std::string &extendpath )
423 return pwrite( name, path, extendpath, type_ );
426 std::string write (
const std::string &name,
427 VTK::OutputType type,
431 addPartitionData( rank );
432 return vtkWriter_->write( name, type, rank,
size );
435 std::string write (
const std::string &name,
439 return write( name, type_, rank,
size );
443 const GridPartType &gridPart_;
444 VTKWriterType *vtkWriter_;
446 VTK::OutputType type_;
454 template<
class Gr
idPart,
bool subsampling >
455 class VTKIOBase< GridPart, subsampling >::VTKWriter
462 using BaseType::write;
463 using BaseType::pwrite;
466 VTKWriter(
const GridPartType &
gridPart,
467 VTK::DataMode dm = VTK::conforming )
477 template<
class Gr
idPart,
bool subsampling >
478 class VTKIOBase< GridPart, subsampling >::SubsamplingVTKWriter
485 using BaseType::write;
486 using BaseType::pwrite;
489 SubsamplingVTKWriter(
const GridPartType &
gridPart,
491 bool coerceToSimplex =
false )
492 : BaseType(
gridPart, intervals, coerceToSimplex )
501 template<
class Gr
idPart >
502 class VTKIO< GridPart, false >
503 :
public VTKIOBase< GridPart, false >
505 typedef VTKIO< GridPart > ThisType;
508 typedef typename BaseType::VTKWriterType VTKWriterType;
511 typedef typename BaseType::GridPartType GridPartType;
513 VTKIO (
const GridPartType &gridPart, VTK::DataMode dm,
const ParameterReader ¶meter = Parameter::container() )
514 : BaseType( gridPart, new VTKWriterType( gridPart, dm ), parameter )
517 explicit VTKIO (
const GridPartType &gridPart,
const ParameterReader ¶meter = Parameter::container() )
518 : VTKIO( gridPart, VTK::
conforming, parameter )
527 template<
class Gr
idPart >
528 class VTKIO< GridPart, true >
529 :
public VTKIOBase< GridPart , true >
531 typedef VTKIO< GridPart, true > ThisType;
534 typedef typename BaseType::VTKWriterType VTKWriterType;
537 typedef typename BaseType::GridPartType GridPartType;
539 explicit VTKIO (
const GridPartType &gridPart,
Dune::RefinementIntervals intervals,
bool coerceToSimplex,
const ParameterReader ¶meter = Parameter::container() )
540 : BaseType( gridPart, new VTKWriterType( gridPart, intervals, coerceToSimplex ), parameter )
543 explicit VTKIO (
const GridPartType &gridPart,
unsigned int level,
bool coerceToSimplex,
const ParameterReader ¶meter = Parameter::container() )
547 VTKIO (
const GridPartType &gridPart,
unsigned int level,
const ParameterReader ¶meter = Parameter::container() )
548 : VTKIO( gridPart, level, false, parameter )
551 [[deprecated(
"pass level as unsigned int" )]]
552 VTKIO (
const GridPartType &gridPart,
int level,
const ParameterReader ¶meter = Parameter::container() )
553 : VTKIO( gridPart, level, false, parameter )
556 VTKIO (
const GridPartType &gridPart,
const ParameterReader ¶meter = Parameter::container() )
557 : VTKIO( gridPart, 0, false, parameter )
560 VTKIO (
const GridPartType &gridPart,
bool coerceToSimplex,
const ParameterReader ¶meter = Parameter::container() )
561 : VTKIO( gridPart, 0, coerceToSimplex, parameter )
570 template<
class Gr
idPart >
571 using SubsamplingVTKIO = VTKIO< GridPart, true >;
consecutive, persistent index set for the leaf level based on the grid's hierarchy index set
Definition: adaptiveleafindexset.hh:1351
BaseType::GridPartType GridPartType
type of the underlying grid part
Definition: discretefunction.hh:609
DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType
type of associated function space
Definition: discretefunction.hh:101
FunctionSpaceTraits::RangeType RangeType
Type of range vector (using type of range field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:71
@ dimDomain
dimension of domain vector space
Definition: functionspaceinterface.hh:46
@ dimRange
dimension of range vector space
Definition: functionspaceinterface.hh:48
FunctionSpaceTraits::DomainType DomainType
Type of domain vector (using type of domain field) has a Dune::FieldVector type interface.
Definition: functionspaceinterface.hh:67
A vector valued function space.
Definition: functionspace.hh:60
forward declaration
Definition: discretefunction.hh:51
Output using VTK.
Definition: vtkio.hh:134
const GridPartType & gridPart() const
return grid part
Definition: vtkio.hh:292
Holds the number of refined intervals per axis needed for virtual and static refinement.
Definition: base.cc:94
Writer for the output of subsampled grid functions in the vtk format.
Definition: subsamplingvtkwriter.hh:40
A base class for grid functions with any return type and dimension.
Definition: function.hh:42
virtual double evaluate(int comp, const Entity &e, const Dune::FieldVector< ctype, dim > &xi) const=0
evaluate single component comp in the entity e at local coordinates xi
Writer for the output of grid functions in the vtk format.
Definition: vtkwriter.hh:95
Definition of the DUNE_NO_DEPRECATED_* macros.
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:43
@ conforming
Output conforming data.
Definition: common.hh:73
RefinementIntervals refinementLevels(int levels)
Creates a RefinementIntervals object.
Definition: base.cc:117
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
Static tag representing a codimension.
Definition: dimension.hh:24
Provides subsampled file i/o for the visualization toolkit.
Provides file i/o for the visualization toolkit.