1#ifndef DUNE_FEM_VTKIO_HH
2#define DUNE_FEM_VTKIO_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>
25 template<
class Gr
idPart,
bool subsampling = false >
35 :
public VTKFunction< typename DF::GridPartType::GridViewType >
46 static const int dimRange = FunctionSpaceType::dimRange;
47 static const int dimDomain = FunctionSpaceType::dimDomain;
49 typedef typename FunctionSpaceType::DomainType
DomainType;
50 typedef typename FunctionSpaceType::RangeType
RangeType;
53 typedef typename GridPartType::template Codim< 0 >::EntityType
EntityType;
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 )
75 return (!vector_) ? 1 : 3;
82 localFunction_.bind( e );
83 typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
85 localFunction_.evaluate(xi,val);
86 localFunction_.unbind();
88 RangeFieldType outVal( 0 );
92 outVal = val[ comp + component_ ] ;
95 outVal = val[ component_ ] ;
103 virtual std::string
name ()
const
105 std::stringstream ret;
112 ret <<
"_vec" << component_;
120 const std::string name_ ;
123 const int component_;
132 template<
class Gr
idPart,
bool subsampling >
143 typedef typename std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
154 :
public VTKFunction< GridViewType >
159 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
166 const std::string&
name,
167 const int rank,
const int nThreads )
168 : iterators_(
gridPart ), name_(
name ), rank_( rank ), nThreads_( nThreads ) {}
174 virtual int ncomps ()
const {
return 1; }
180 const int thread = iterators_.
thread( e );
181 return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
185 virtual std::string
name ()
const
192 const std::string name_;
199 :
public VTKFunction< GridViewType >
204 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
216 virtual int ncomps ()
const {
return 1; }
222 return e.geometry().volume();
226 virtual std::string
name ()
const
228 return std::string(
"volume");
235 const std::string names[] = {
"none",
"rank",
"rank+thread",
"rank/thread" };
236 return parameter.getEnum(
"fem.io.partitioning", names, 0 );
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 ) ];
254 std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
257 const int rank = ( myRank < 0 ) ?
gridPart_.comm().rank() : myRank ;
261 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(
gridPart_,
"rank", rank, nThreads) );
267 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(
gridPart_,
"rank", rank, -1) );
270 std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(
gridPart_,
"thread", 0, nThreads) );
277 template <
class DF >
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;
300 static const int dimRange = DF::FunctionSpaceType::dimRange;
301 for(
int i = 0;i < dimRange; ++i )
303 if ( notComplex<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,
334 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
337 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
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>() )
370 const std::string& dataName =
"" ,
373 if ( notComplex<DF>() )
375 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
381 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
384 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
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 );
405 std::string
write (
const std::string &name )
410 std::string
pwrite (
const std::string &name,
411 const std::string &
path,
412 const std::string &extendpath,
413 VTK::OutputType type )
419 std::string
pwrite (
const std::string &name,
420 const std::string &
path,
421 const std::string &extendpath )
426 std::string
write (
const std::string &name,
427 VTK::OutputType type,
432 return vtkWriter_->write( name, type, rank, size );
435 std::string
write (
const std::string &name,
454 template<
class Gr
idPart,
bool subsampling >
456 :
public Dune::VTKWriter< GridViewType >
458 typedef Dune::VTKWriter< GridViewType > BaseType;
462 using BaseType::write;
463 using BaseType::pwrite;
467 VTK::DataMode dm = VTK::conforming )
468 : BaseType(
gridPart.gridView(), dm )
477 template<
class Gr
idPart,
bool subsampling >
479 :
public Dune::SubsamplingVTKWriter< GridViewType >
481 typedef Dune::SubsamplingVTKWriter< GridViewType > BaseType;
485 using BaseType::write;
486 using BaseType::pwrite;
490 Dune::RefinementIntervals intervals,
491 bool coerceToSimplex =
false )
492 : BaseType(
gridPart.gridView(), intervals, coerceToSimplex )
501 template<
class Gr
idPart >
514 :
BaseType( gridPart, new VTKWriterType( gridPart, dm ), parameter )
518 :
VTKIO( gridPart, VTK::conforming, parameter )
527 template<
class Gr
idPart >
540 :
BaseType( gridPart, new VTKWriterType( gridPart, intervals, coerceToSimplex ), parameter )
544 :
VTKIO( gridPart,
Dune::refinementLevels(level), coerceToSimplex, parameter )
548 :
VTKIO( gridPart, level, false, parameter )
552 :
VTKIO( gridPart, level, false, parameter )
556 : VTKIO( gridPart, 0, false, parameter )
560 : VTKIO( gridPart, 0, coerceToSimplex, parameter )
569 template<
class Gr
idPart >
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
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 ¶meter=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 ¶meter=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
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
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
VTKIO(const GridPartType &gridPart, Dune::RefinementIntervals intervals, bool coerceToSimplex, const ParameterReader ¶meter=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