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>
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 )
62 : localFunction_( df ),
63 name_( ( dataName.size() > 0 ) ? dataName : df.
name() ),
65 typeOfField_( typeOfField ),
66 component_( component )
76 return (!vector_) ? 1 : 3;
83 localFunction_.bind( e );
84 typedef typename LocalFunctionType::RangeFieldType RangeFieldType;
86 localFunction_.evaluate(xi,val);
87 localFunction_.unbind();
89 RangeFieldType outVal( 0 );
93 outVal = val[ comp + component_ ] ;
96 outVal = val[ component_ ] ;
97 if (typeOfField_ ==
TypeOfField::real || typeOfField_ == TypeOfField::complex_real )
104 virtual std::string
name ()
const
106 std::stringstream ret;
108 if (typeOfField_ == TypeOfField::complex_real)
110 if (typeOfField_ == TypeOfField::complex_imag)
113 ret <<
"_vec" << component_;
121 const std::string name_ ;
124 const int component_;
133 template<
class Gr
idPart,
bool subsampling >
144 typedef typename std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type
155 :
public VTKFunction< GridViewType >
160 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
167 const std::string&
name,
168 const int rank,
const int nThreads )
169 : iterators_(
gridPart ), name_(
name ), rank_( rank ), nThreads_( nThreads ) {}
175 virtual int ncomps ()
const {
return 1; }
181 const int thread = iterators_.
thread( e );
182 return (nThreads_ < 0) ? double( rank_ ) : double( rank_ * nThreads_ + thread );
186 virtual std::string
name ()
const
193 const std::string name_;
200 :
public VTKFunction< GridViewType >
205 typedef typename GridViewType :: template Codim< 0 >::Entity
EntityType;
217 virtual int ncomps ()
const {
return 1; }
223 return e.geometry().volume();
227 virtual std::string
name ()
const
229 return std::string(
"volume");
236 const std::string names[] = {
"none",
"rank",
"rank+thread",
"rank/thread" };
237 return parameter.getEnum(
"fem.io.partitioning", names, 0 );
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 ) ];
255 std::shared_ptr<VolumeData> volumePtr( std::make_shared<VolumeData>() );
258 const int rank = ( myRank < 0 ) ?
gridPart_.comm().rank() : myRank ;
262 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(
gridPart_,
"rank", rank, nThreads) );
268 std::shared_ptr<PartitioningData> dataRankPtr( std::make_shared<PartitioningData>(
gridPart_,
"rank", rank, -1) );
271 std::shared_ptr<PartitioningData> dataThreadPtr( std::make_shared<PartitioningData>(
gridPart_,
"thread", 0, nThreads) );
278 template <
class DF >
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;
301 static const int dimRange = DF::FunctionSpaceType::dimRange;
302 for(
int i = 0;i < dimRange; ++i )
304 if ( notComplex<DF>() )
324 const std::string& dataName =
"" ,
327 if ( notComplex<DF>() )
329 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
335 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
338 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
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 )
351 if ( notComplex<DF>() )
371 const std::string& dataName =
"" ,
374 if ( notComplex<DF>() )
376 std::shared_ptr<VTKFunctionWrapper< DF > > ptr( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
382 std::shared_ptr<VTKFunctionWrapper< DF > > ptrR( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
385 std::shared_ptr<VTKFunctionWrapper< DF > > ptrI( std::make_shared<
VTKFunctionWrapper< DF > >( df, dataName, startPoint,
396 std::string
write (
const std::string &name, VTK::OutputType type )
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 );
406 std::string
write (
const std::string &name )
411 std::string
pwrite (
const std::string &name,
412 const std::string &
path,
413 const std::string &extendpath,
414 VTK::OutputType type )
420 std::string
pwrite (
const std::string &name,
421 const std::string &
path,
422 const std::string &extendpath )
427 std::string
write (
const std::string &name,
428 VTK::OutputType type,
433 return vtkWriter_->write( name, type, rank, size );
436 std::string
write (
const std::string &name,
455 template<
class Gr
idPart,
bool subsampling >
457 :
public Dune::VTKWriter< GridViewType >
459 typedef Dune::VTKWriter< GridViewType > BaseType;
463 using BaseType::write;
464 using BaseType::pwrite;
468 VTK::DataMode dm = VTK::conforming )
469 : BaseType(
gridPart.gridView(), dm )
478 template<
class Gr
idPart,
bool subsampling >
480 :
public Dune::SubsamplingVTKWriter< GridViewType >
482 typedef Dune::SubsamplingVTKWriter< GridViewType > BaseType;
486 using BaseType::write;
487 using BaseType::pwrite;
491 Dune::RefinementIntervals intervals,
492 bool coerceToSimplex =
false )
493 : BaseType(
gridPart.gridView(), intervals, coerceToSimplex )
502 template<
class Gr
idPart >
515 :
BaseType( gridPart, new VTKWriterType( gridPart, dm ), parameter )
519 :
VTKIO( gridPart, VTK::conforming, parameter )
528 template<
class Gr
idPart >
541 :
BaseType( gridPart, new VTKWriterType( gridPart, intervals, coerceToSimplex ), parameter )
545 :
VTKIO( gridPart,
Dune::refinementLevels(level), coerceToSimplex, parameter )
549 :
VTKIO( gridPart, level, false, parameter )
553 :
VTKIO( gridPart, level, false, parameter )
557 : VTKIO( gridPart, 0, false, parameter )
561 : VTKIO( gridPart, 0, coerceToSimplex, parameter )
570 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
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
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
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
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: vtkio.hh:54
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 ¶meter=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
std::conditional< subsampling, SubsamplingVTKWriter, VTKWriter >::type VTKWriterType
Definition: vtkio.hh:142
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
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
const GridPartType & gridPart() const
return grid part
Definition: vtkio.hh:293
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 ¶meter=Parameter::container()) const
Definition: vtkio.hh:233
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
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:160
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
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
GridViewType::template Codim< 0 >::Entity EntityType
Definition: vtkio.hh:205
VolumeData()
constructor taking discrete function
Definition: vtkio.hh:211
virtual std::string name() const
get name
Definition: vtkio.hh:227
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
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
VTKIO(const GridPartType &gridPart, Dune::RefinementIntervals intervals, bool coerceToSimplex, const ParameterReader ¶meter=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