1 #ifndef DUNE_FEM_VTKIO_HH 2 #define DUNE_FEM_VTKIO_HH 4 #include <dune/common/typetraits.hh> 6 #include <dune/common/deprecated.hh> 8 #include <dune/grid/io/file/vtk/vtkwriter.hh> 9 #include <dune/grid/io/file/vtk/subsamplingvtkwriter.hh> 26 template<
class Gr
idPart,
bool subsampling = false >
36 :
public VTKFunction< typename DF::GridPartType::GridViewType >
47 static const int dimRange = FunctionSpaceType::dimRange;
48 static const int dimDomain = FunctionSpaceType::dimDomain;
50 typedef typename FunctionSpaceType::DomainType
DomainType;
51 typedef typename FunctionSpaceType::RangeType
RangeType;
54 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
60 const std::string& dataName,
61 int component,
bool vector,
TypeOfField typeOfField )
63 name_( ( dataName.size() > 0 ) ? dataName : df.
name() ),
65 typeOfField_( typeOfField ),
66 component_( component )
76 return (!vector_) ? 1 : 3;
81 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const 83 typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
84 const LocalFunctionType lf = discFunc_.localFunction(e);
87 RangeFieldType outVal( 0 );
90 if( comp <= dimDomain )
91 outVal = val[ comp + component_ ] ;
94 outVal = val[ component_ ] ;
95 if (typeOfField_ ==
TypeOfField::real || typeOfField_ == TypeOfField::complex_real )
102 virtual std::string
name ()
const 104 std::stringstream ret;
106 if (typeOfField_ == TypeOfField::complex_real)
108 if (typeOfField_ == TypeOfField::complex_imag)
111 ret <<
"_vec" << component_;
118 const DiscreteFunctionType &discFunc_;
119 const std::string name_ ;
122 const int component_;
131 template<
class Gr
idPart,
bool subsampling >
140 class SubsamplingVTKWriter;
142 typedef typename conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
153 :
public VTKFunction< GridViewType >
158 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
165 const std::string&
name,
166 const int rank,
const int nThreads )
167 : iterators_( gridPart ), name_( name ), rank_( rank ), nThreads_( nThreads ) {}
173 virtual int ncomps ()
const {
return 1; }
177 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const 179 const int thread = iterators_.thread( e );
180 return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
184 virtual std::string
name ()
const 190 ThreadIteratorType iterators_;
191 const std::string name_;
198 :
public VTKFunction< GridViewType >
203 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
215 virtual int ncomps ()
const {
return 1; }
219 virtual double evaluate (
int comp,
const EntityType &e,
const LocalCoordinateType &xi )
const 221 return e.geometry().volume();
225 virtual std::string
name ()
const 227 return std::string(
"volume");
234 const std::string names[] = {
"none",
"rank",
"rank+thread",
"rank/thread" };
235 return parameter.getEnum(
"fem.io.partitioning", names, 0 );
240 : gridPart_( gridPart ),
241 vtkWriter_( vtkWriter ),
242 addPartition_( getPartitionParameter( parameter ) )
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 ) ];
251 if( addPartition_ > 0 )
253 std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
254 vtkWriter_->addCellData( volumePtr );
256 const int rank = ( myRank < 0 ) ? gridPart_.comm().rank() : myRank ;
258 if( addPartition_ <= 2 )
260 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_,
"rank", rank, nThreads) );
261 vtkWriter_->addCellData( dataRankPtr );
266 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(gridPart_,
"rank", rank, -1) );
267 vtkWriter_->addCellData( dataRankPtr );
269 std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(gridPart_,
"thread", 0, nThreads) );
270 vtkWriter_->addCellData( dataThreadPtr );
276 template <
class DF >
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;
299 static const int dimRange = DF::FunctionSpaceType::dimRange;
302 if ( notComplex<DF>() )
306 vtkWriter_->addCellData( ptr );
312 vtkWriter_->addCellData( ptrR );
315 vtkWriter_->addCellData( ptrI );
322 const std::string& dataName =
"" ,
325 if ( notComplex<DF>() )
327 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
329 vtkWriter_->addCellData( ptr );
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 );
345 static const int dimRange = DF::FunctionSpaceType::dimRange;
346 std::string
name = ( dataName.size() > 0 ) ? dataName : df.name() ;
349 if ( notComplex<DF>() )
353 vtkWriter_->addVertexData( ptr );
359 vtkWriter_->addVertexData( ptrR );
362 vtkWriter_->addVertexData( ptrI );
369 const std::string& dataName =
"" ,
372 if ( notComplex<DF>() )
374 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
376 vtkWriter_->addVertexData( ptr );
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 );
394 std::string
write (
const std::string &
name, VTK::OutputType type )
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 );
401 return vtkWriter_->write( name, type );
406 return write( name, type_ );
410 const std::string &
path,
411 const std::string &extendpath,
412 VTK::OutputType type )
415 return vtkWriter_->pwrite( name, path, extendpath, type );
419 const std::string &
path,
420 const std::string &extendpath )
422 return pwrite( name, path, extendpath, type_ );
426 VTK::OutputType type,
430 addPartitionData( rank );
431 return vtkWriter_->write( name, type, rank, size );
438 return write( name, type_, rank, size );
453 template<
class Gr
idPart,
bool subsampling >
455 :
public Dune::VTKWriter< GridViewType >
457 typedef Dune::VTKWriter< GridViewType > BaseType;
461 using BaseType::write;
462 using BaseType::pwrite;
466 VTK::DataMode dm = VTK::conforming )
467 : BaseType( static_cast< GridViewType > ( gridPart ), dm )
476 template<
class Gr
idPart,
bool subsampling >
478 :
public Dune::SubsamplingVTKWriter< GridViewType >
480 typedef Dune::SubsamplingVTKWriter< GridViewType > BaseType;
484 using BaseType::write;
485 using BaseType::pwrite;
490 bool coerceToSimplex =
false )
491 : BaseType( static_cast< GridViewType > ( gridPart ),
492 level, coerceToSimplex )
501 template<
class Gr
idPart >
508 typedef typename BaseType::VTKWriterType VTKWriterType;
514 : BaseType( gridPart, new VTKWriterType( gridPart, dm ), parameter )
518 :
VTKIO( gridPart, VTK::conforming, parameter )
527 template<
class Gr
idPart >
534 typedef typename BaseType::VTKWriterType VTKWriterType;
540 : BaseType( gridPart, new VTKWriterType( gridPart, level, coerceToSimplex ), parameter )
544 :
VTKIO( gridPart, level, false, parameter )
548 :
VTKIO( gridPart, level,
false, parameter )
552 :
VTKIO( gridPart, 0,
false, parameter )
556 :
VTKIO( gridPart, 0, coerceToSimplex, parameter )
566 template<
class Gr
idPart >
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 ¶meter=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 ¶meter=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
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
virtual double evaluate(int comp, const EntityType &e, const LocalCoordinateType &xi) const
Definition: vtkio.hh:219
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 ¶meter=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
const GridPartType & gridPart() const
return grid part
Definition: vtkio.hh:291
DiscreteFunctionType::LocalFunctionType LocalFunctionType
Definition: vtkio.hh:44
/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