1 #ifndef DUNE_FEM_DATAOUTPUT_HH 2 #define DUNE_FEM_DATAOUTPUT_HH 5 #define USE_VTKWRITER 1 22 #ifndef ENABLE_VTXPROJECTION 23 #define ENABLE_VTXPROJECTION 1 28 #if ENABLE_VTXPROJECTION 31 #endif // #if ENABLE_VTXPROJECTION 32 #endif // #if USE_VTKWRITER 36 #define USE_GRAPE HAVE_GRAPE 40 #include <dune/grid/io/visual/grapedatadisplay.hh> 59 :
public LocalParameter< DataOutputParameters, DataOutputParameters >
72 : keyPrefix_(
"fem.io." ), parameter_(
parameter )
76 virtual std::string
path()
const 83 std::string absPath =
parameter().
getValue< std::string >(
"fem.prefix",
"." ) +
"/";
84 const std::string relPath =
path();
99 static const std::string formatTable[]
100 = {
"vtk-cell",
"vtk-vertex",
"sub-vtk-cell",
"binary" ,
"gnuplot" ,
"none" };
101 int format =
parameter().
getEnum( keyPrefix_ +
"outputformat", formatTable, 1 );
114 return (
parameter().getValue( keyPrefix_ +
"grapedisplay", 0 ) == 1);
182 template<
class Gr
idImp,
class DataImp >
189 template< class Grid, class OutputTuple, int N = std::tuple_size< OutputTuple >::value >
192 template<
class Gr
id,
class OutputTuple,
int N >
195 typedef typename TypeTraits< typename std::tuple_element< 0, OutputTuple >::type >::PointeeType
DFType;
196 typedef typename DFType :: DiscreteFunctionSpaceType :: GridPartType
GridPartType;
199 : gridPart_( getGridPart( data ) )
210 const DFType *df = std::get< 0 >( data );
212 return df->space().gridPart();
219 template<
class Gr
id,
class OutputTuple >
225 : gridPart_( const_cast< Grid & >( grid ) )
239 template<
class VTKIOType >
244 virtual void add ( VTKIOType & )
const = 0;
251 #if ENABLE_VTXPROJECTION 252 template<
class VTKIOType,
class DFType >
255 template<
class VTKOut >
257 #endif // #if ENABLE_VTXPROJECTION 259 template<
class VTKOut >
261 #endif // #if USE_VTKWRITER 263 template<
class Gr
idPartType >
267 enum OutputFormat { vtk = 0, vtkvtx = 1, subvtk = 2 , binary = 3, gnuplot = 4, none = 5 };
309 pvd_ <<
" </Collection>" << std::endl;
310 pvd_ <<
"</VTKFile>" << std::endl;
325 consistentSaveStep( tp );
327 const double saveTimeEps = 2*std::numeric_limits< double >::epsilon()*saveTime_;
328 const bool writeStep = (saveStep_ > 0) && (tp.
time() - saveTime_ >= -saveTimeEps);
329 const bool writeCount = (saveCount_ > 0) && (writeCalls_ % saveCount_ == 0);
330 return param_->willWrite( writeStep || writeCount );
336 return param_->willWrite( (saveCount_ > 0) && (writeCalls_ % saveCount_ == 0) );
342 void write(
const std::string& outstring )
const 345 writeData( writeCalls_, outstring );
353 writeData( writeCalls_,
"" );
364 writeData( tp.
time(), outstring );
380 void writeData (
double sequenceStamp,
const std::string &outstring )
const;
387 writeData( sequenceStamp,
"" );
397 const std::string &
path ()
const 422 std::string writeVTKOutput ()
const;
425 std::string writeGnuPlotOutput ()
const;
430 DUNE_THROW(NotImplemented,
"Format 'binary' not supported by DataOutput, use DataWriter instead." );
437 grapeDisplay( data_ );
442 template<
class OutputTupleType >
443 void grapeDisplay ( OutputTupleType &data )
const;
483 #if ENABLE_VTXPROJECTION 484 template<
class Gr
idImp,
class DataImp >
485 template<
class VTKIOType,
class DFType >
489 typedef typename VTKIOType::GridPartType GridPartType;
490 typedef typename DFType::DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType;
496 VTKFunc (
const GridPartType &gridPart,
const DFType &df )
498 space_( const_cast< GridPartType & >( gridPart ) ),
509 virtual void add ( VTKIOType &vtkio )
const 511 func_ =
new NewFunctionType( df_.name()+
"vtx-prj" , space_ );
512 if( df_.space().continuous() )
522 vtkio.addVertexData( *func_ );
529 LagrangeSpaceType space_;
530 mutable NewFunctionType *func_;
532 #endif // #if ENABLE_VTXPROJECTION 533 #endif // #if USE_VTKWRITER 541 template<
class Gr
idImp,
class DataImp >
542 template<
class VTKOut >
548 conforming_( conforming )
552 template<
class DFType >
557 if( conforming_ || (df->space().order() == 0) )
558 vtkOut_.addCellData( *df );
560 vtkOut_.addVertexData( *df );
564 template<
class OutputTuple >
567 ForEachValue< OutputTuple > forEach( data );
568 forEach.apply( *
this );
575 #endif // #if USE_VTKWRITER 583 #if ENABLE_VTXPROJECTION 584 template<
class Gr
idImp,
class DataImp >
585 template<
class VTKOut >
597 for(
auto& entry : vec_ )
605 template<
class DFType >
611 entry->
add( vtkOut_ );
612 vec_.push_back( entry );
616 template<
class OutputTuple >
619 ForEachValue< OutputTuple > forEach( data );
620 forEach.apply( *
this );
626 std::vector< VTKListEntryType * > vec_;
628 #endif // #if ENABLE_VTXPROJECTION 629 #endif // #if USE_VTKWRITER 636 template<
class Gr
idImp,
class DataImp >
637 template<
class Gr
idPartType >
640 typedef typename GridPartType::template Codim< 0 >::IteratorType::Entity Entity;
652 : out_(out), quad_(quad), i_(i), en_(en)
656 template <
class DFType>
661 auto lf = df->localFunction(en_);
662 typename DFType::DiscreteFunctionSpaceType::RangeType u;
663 lf.evaluate( quad_[ i_ ], u );
665 DUNE_CONSTEXPR
int dimRange = DFType::DiscreteFunctionSpaceType::dimRange;
666 for(
auto k = 0; k < dimRange; ++k )
667 out_ <<
" " << u[ k ];
671 template<
class OutputTuple >
674 ForEachValue< OutputTuple > forEach( data );
675 forEach.apply( *
this );
684 template<
class Gr
idImp,
class DataImp >
695 outputFormat_(vtkvtx),
696 conformingOutput_( false ),
697 param_(parameter.clone())
704 template<
class Gr
idImp,
class DataImp >
716 outputFormat_(vtkvtx),
717 conformingOutput_( false ),
718 param_(parameter.clone())
724 consistentSaveStep( tp );
727 template<
class Gr
idImp,
class DataImp >
734 const auto oldTime = tp.
time() - saveStep_;
735 while( saveTime_ <= oldTime )
738 saveTime_ += saveStep_;
744 template<
class Gr
idImp,
class DataImp >
756 datapref_ += parameter.
prefix();
759 switch( outputFormat )
761 case 0: outputFormat_ = vtk;
break;
762 case 1: outputFormat_ = vtkvtx;
break;
763 case 2: outputFormat_ = subvtk;
break;
764 case 3: outputFormat_ = binary;
break;
765 case 4: outputFormat_ = gnuplot;
break;
766 case 5: outputFormat_ = none;
break;
768 DUNE_THROW(NotImplemented,
"DataOutput::init: wrong output format");
790 std::string name = path_ +
"/" + datapref_;
792 std::cout <<
"opening file: " << name << std::endl;
793 sequence_.open(name.c_str());
795 std::cout <<
"could not write sequence file" << std::endl;
799 std::string name = path_ +
"/" + datapref_;
801 std::cout <<
"opening file: " << name << std::endl;
802 pvd_.open(name.c_str());
804 std::cout <<
"could not write sequence file" << std::endl;
807 pvd_ <<
"<?xml version=\"1.0\"?>" << std::endl;
808 pvd_ <<
"<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
809 pvd_ <<
" <Collection>" << std::endl;
820 template<
class Gr
idImp,
class DataImp >
822 ::writeData (
double sequenceStamp,
const std::string &outstring )
const 824 std::string filename;
827 switch( outputFormat_ )
831 case binary: writeBinaryData( sequenceStamp );
break;
837 filename = writeVTKOutput();
840 DUNE_THROW(NotImplemented,
"DataOutput::write: VTKWriter was disabled by USE_VTKWRITER 0");
841 #endif // #if USE_VTKWRITER 842 case gnuplot : filename = writeGnuPlotOutput();
break;
844 DUNE_THROW(NotImplemented,
"DataOutput::write: wrong output format = " << outputFormat_);
847 if( outputFormat_ != none )
850 sequence_ << writeStep_ <<
" " 861 std::string basefilename;
862 auto pos = filename.find_last_of(
'/' );
863 if (pos == filename.npos)
865 basefilename = filename.substr( pos+1, filename.npos );
866 pvd_ <<
" <DataSet timestep=\"" << sequenceStamp <<
"\" " 867 <<
"group=\"\" part=\"0\" " 868 <<
"file=\""<<basefilename<<
"\"/>" << std::endl;
875 std::cout << myClassName() <<
"[" << grid_.comm().rank() <<
"]::write data" 876 <<
" writestep=" << writeStep_
877 <<
" sequenceStamp=" << sequenceStamp
883 saveTime_ += saveStep_;
889 template<
class Gr
idImp,
class DataImp >
892 std::string filename;
894 const bool vertexData = (outputFormat_ == vtkvtx);
897 const bool parallel = (grid_.comm().size() > 1);
900 auto name =
generateFilename( (parallel ? datapref_ : path_ +
"/" + datapref_ ), writeStep_ );
904 #if ENABLE_VTXPROJECTION 907 GridPartGetterType gp( grid_, data_ );
908 const auto &gridPart = gp.gridPart();
912 VTKIOType vtkio ( gridPart, VTK::conforming, param_->parameter() );
921 filename = vtkio.pwrite( name, path_,
"." );
926 filename = vtkio.write( name );
930 else if ( outputFormat_ == vtk )
933 GridPartGetterType gp( grid_, data_ );
937 VTKIOType vtkio( gp.gridPart(), conformingOutput_ ? VTK::conforming : VTK::nonconforming, param_->parameter() );
947 filename = vtkio.pwrite( name, path_,
"." );
952 filename = vtkio.write( name );
955 else if ( outputFormat_ == subvtk )
958 GridPartGetterType gp( grid_, data_ );
962 VTKIOType vtkio ( gp.gridPart(),
static_cast< unsigned int >( param_->subsamplingLevel() ), param_->parameter() );
972 filename = vtkio.pwrite( name, path_,
"." );
977 filename = vtkio.write( name );
982 #endif // #if USE_VTKWRITER 985 template<
class Gr
idImp,
class DataImp >
991 std::ofstream gnuout(name.c_str());
992 gnuout << std::scientific << std::setprecision( 16 );
995 GridPartGetterType gp( grid_, data_ );
996 const auto &gridPart = gp.gridPart();
999 typedef typename GridPartGetterType::GridPartType GridPartType;
1000 typedef typename GridPartType::GridViewType GridViewType;
1001 for(
const auto entity : elements( static_cast<GridViewType>( gridPart ) ) )
1004 for(
unsigned int i = 0; i < quad.nop(); ++i )
1006 const auto x = entity.geometry().global( quad.point( i ) );
1007 for(
auto k = 0; k < x.dimension; ++k )
1008 gnuout << (k > 0 ?
" " :
"") << x[ k ];
1011 gnuout << std::endl;
1019 template<
class Gr
idImp,
class DataImp >
1020 template<
class OutputTupleType >
1023 GrapeDataDisplay<GridType> grape(grid_);
1027 #else // #if USE_GRAPE 1028 template<
class Gr
idImp,
class DataImp >
1029 template<
class OutputTupleType >
1032 std::cerr <<
"WARNING: HAVE_GRAPE == 0, but grapeDisplay == true, recompile with grape support for online display!" << std::endl;
1034 #endif // #else // #if USE_GRAPE 1040 #endif // #ifndef DUNE_FEM_DATAOUTPUT_HH DataOutput(const GridType &grid, OutPutDataType &data, const ParameterReader ¶meter=Parameter::container())
Definition: dataoutput.hh:283
static void project(const Function &f, DiscreteFunction &u, Weight &weight)
Definition: vtxprojection.hh:32
std::string writeGnuPlotOutput() const
Definition: dataoutput.hh:986
void visit(DFType *df)
Applies the setting on every DiscreteFunction/LocalFunction pair.
Definition: dataoutput.hh:606
virtual int subsamplingLevel() const
save data every subsamplingLevel (fem.io.subsamplinglevel)
Definition: dataoutput.hh:133
std::string path_
Definition: dataoutput.hh:450
DataImp OutPutDataType
type of data tuple
Definition: dataoutput.hh:272
void visit(DFType *df)
Applies the setting on every DiscreteFunction/LocalFunction pair.
Definition: dataoutput.hh:553
VTKOutputerDG(VTKOut &vtkOut, bool conforming=false)
Constructor.
Definition: dataoutput.hh:546
virtual bool willWrite(const TimeProviderBase &tp) const
returns true if data will be written on next write call
Definition: dataoutput.hh:322
OutputFormat outputFormat_
Definition: dataoutput.hh:466
DFType::DiscreteFunctionSpaceType::GridPartType GridPartType
Definition: dataoutput.hh:196
GnuplotOutputer(std::ostream &out, CachingQuadrature< GridPartType, 0 > &quad, int i, const Entity &en)
Constructor.
Definition: dataoutput.hh:648
std::string writeVTKOutput() const
Definition: dataoutput.hh:890
int writeCalls() const
return write calls
Definition: dataoutput.hh:409
std::string datapref_
Definition: dataoutput.hh:451
const ParameterReader & parameter() const noexcept
Definition: dataoutput.hh:170
virtual ~VTKListEntry()
Definition: dataoutput.hh:242
double saveStep_
Definition: dataoutput.hh:462
void forEach(OutputTuple &data)
Definition: dataoutput.hh:672
virtual double startsavetime() const
value of first save time (no parameter available)
Definition: dataoutput.hh:151
int writeStep_
Definition: dataoutput.hh:456
virtual void writeBinaryData(const double) const
write binary data
Definition: dataoutput.hh:428
int getEnum(const std::string &key, const std::string(&values)[n]) const
Definition: reader.hh:215
int writeStep() const
return write step
Definition: dataoutput.hh:403
void write(const TimeProviderBase &tp) const
write given data to disc, evaluates parameter savecount and savestep
Definition: dataoutput.hh:371
Implementation of the Dune::Fem::IOInterface. This class manages data output. Available output format...
Definition: dataoutput.hh:183
~VTKOutputerLagrange()
destructor, delete entries of list
Definition: dataoutput.hh:595
virtual bool willWrite(bool write) const
Definition: dataoutput.hh:158
Definition: adaptivefunction/adaptivefunction.hh:34
virtual std::string absolutePath() const
Definition: dataoutput.hh:81
DataOutputParameters(std::string keyPrefix, const ParameterReader ¶meter=Parameter::container())
Definition: dataoutput.hh:67
DataOutput(const GridType &grid, OutPutDataType &data, const TimeProviderBase &tp, const ParameterReader ¶meter=Parameter::container())
Definition: dataoutput.hh:296
virtual int startcall() const
number of first call (no parameter available)
Definition: dataoutput.hh:145
void write(const TimeProviderBase &tp, const std::string &outstring) const
write given data to disc, evaluates parameter savecount and savestep
Definition: dataoutput.hh:361
void visit(DFType *df)
Applies the setting on every DiscreteFunction/LocalFunction pair.
Definition: dataoutput.hh:657
const GridType & grid_
type of this class
Definition: dataoutput.hh:446
Definition: dataoutput.hh:260
virtual std::string prefix() const
base of file name for data file (fem.io.datafileprefix)
Definition: dataoutput.hh:91
void write(const std::string &outstring) const
write given data to disc, evaluates parameter savecount
Definition: dataoutput.hh:342
const DataOutputParameters * param_
Definition: dataoutput.hh:470
bool conformingOutput_
Definition: dataoutput.hh:467
OutputFormat
Definition: dataoutput.hh:267
virtual bool conformingoutput() const
Definition: dataoutput.hh:105
const GridPartType & gridPart() const
Definition: dataoutput.hh:202
Definition: vtxprojection.hh:17
Parameter class for Dune::Fem::DataOutput.
Definition: dataoutput.hh:57
virtual int startcounter() const
number for first data file (no parameter available)
Definition: dataoutput.hh:139
~VTKFunc()
Definition: dataoutput.hh:504
Definition: dataoutput.hh:253
general base for time providers
Definition: timeprovider.hh:35
Definition: dataoutput.hh:190
std::string generateFilename(const std::string &fn, int ntime, int precision=6)
Definition: iointerface.hh:47
Definition: dataoutput.hh:264
VTKListEntry()
Definition: dataoutput.hh:247
virtual void display() const
display data with grape
Definition: dataoutput.hh:434
virtual bool grapedisplay() const
use online grape display (fem.io.grapedisplay)
Definition: dataoutput.hh:111
Definition: io/parameter.hh:549
Definition: coordinate.hh:4
DataOutputParameters(const ParameterReader ¶meter=Parameter::container())
Definition: dataoutput.hh:71
int writeCalls_
Definition: dataoutput.hh:458
Definition: dataoutput.hh:256
const GridPartType & gridPart_
Definition: dataoutput.hh:215
GridImp GridType
type of grid used
Definition: dataoutput.hh:270
Definition: dataoutput.hh:240
virtual ~DataOutput()
destructor
Definition: dataoutput.hh:303
void init(const DataOutputParameters ¶meter)
initialize data writer
Definition: dataoutput.hh:746
static ParameterContainer & container()
Definition: io/parameter.hh:190
AdaptiveLeafGridPart< Grid > GridPartType
Definition: dataoutput.hh:222
static void write(const std::string &filename, const std::string &fileextension="", bool writeAll=true)
write the parameter database to a file
Definition: io/parameter.hh:515
int saveCount_
Definition: dataoutput.hh:464
virtual bool willWrite() const
returns true if data will be written on next write call
Definition: dataoutput.hh:334
virtual void add(VTKIOType &vtkio) const
Definition: dataoutput.hh:509
void consistentSaveStep(const TimeProviderBase &tp) const
Definition: dataoutput.hh:729
OutPutDataType & data_
Definition: dataoutput.hh:447
GridPartGetter(const Grid &, const OutputTuple &data)
Definition: dataoutput.hh:198
void writeData(double sequenceStamp) const
write data with a given sequence stamp
Definition: dataoutput.hh:385
static void createGlobalPath(const CommunicatorType &comm, const std::string &path)
create global path for data output
Definition: iointerface.hh:267
static const GridPartType & getGridPart(const OutputTuple &data)
Definition: dataoutput.hh:208
static bool verbose()
obtain the cached value for fem.verbose
Definition: io/parameter.hh:444
IOInterface to write data to hard disk.
Definition: iointerface.hh:161
void move(ArrayInterface< T > &array, const unsigned int oldOffset, const unsigned int newOffset, const unsigned int length)
Definition: array_inline.hh:38
void forEach(OutputTuple &data)
Definition: dataoutput.hh:565
void write() const
write given data to disc, evaluates parameter savecount
Definition: dataoutput.hh:350
virtual std::string path() const
path where the data is stored (always relative to fem.prefix)
Definition: dataoutput.hh:76
ParameterReader parameter_
Definition: dataoutput.hh:64
const std::string & path() const
return output path name
Definition: dataoutput.hh:397
Definition: adaptiveleafgridpart.hh:248
T getValue(const std::string &key) const
get mandatory parameter
Definition: reader.hh:149
const std::string keyPrefix_
Definition: dataoutput.hh:63
static void addToDisplay(Disp &disp, const DINFO *dinf, double time, Tuple &tuple)
Definition: iotuple.hh:270
static void interpolateFunction(const Function &function, DiscreteFunctionType &discreteFunction)
Definition: lagrangeinterpolation.hh:82
const GridPartType gridPart_
Definition: dataoutput.hh:234
std::ofstream pvd_
Definition: dataoutput.hh:469
TypeTraits< typename std::tuple_element< 0, OutputTuple >::type >::PointeeType DFType
Definition: dataoutput.hh:195
void forEach(OutputTuple &data)
Definition: dataoutput.hh:617
double time() const
obtain the current time
Definition: timeprovider.hh:94
std::ofstream sequence_
Definition: dataoutput.hh:468
GridPartGetter(const Grid &grid, const OutputTuple &)
Definition: dataoutput.hh:224
VTKOutputerLagrange(VTKOut &vtkOut)
Constructor.
Definition: dataoutput.hh:589
VTKFunc(const GridPartType &gridPart, const DFType &df)
Definition: dataoutput.hh:496
double saveTime() const
return save time
Definition: dataoutput.hh:415
virtual int outputformat() const
format of output (fem.io.outputformat)
Definition: dataoutput.hh:97
bool grapeDisplay_
Definition: dataoutput.hh:453
DataOutput(const GridType &grid, OutPutDataType &data, const DataOutputParameters ¶meter)
Constructor creating data output class.
Definition: dataoutput.hh:686
double saveTime_
Definition: dataoutput.hh:460
const GridPartType & gridPart() const
Definition: dataoutput.hh:228
void writeData(double sequenceStamp, const std::string &outstring) const
write data with a given sequence stamp and outstring
Definition: dataoutput.hh:822
void grapeDisplay(OutputTupleType &data) const
display data with grape
Definition: dataoutput.hh:1030
virtual const char * myClassName() const
print class name
Definition: dataoutput.hh:391
virtual double savestep() const
save data every savestep interval (fem.io.savestep)
Definition: dataoutput.hh:121
virtual bool writeMode() const
Definition: dataoutput.hh:165
virtual int savecount() const
save data every savecount calls to write method (fem.io.savecount)
Definition: dataoutput.hh:127