Commit 932d92ee authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

cleanup of reader and writer

parent 44d12215
...@@ -4,23 +4,67 @@ ...@@ -4,23 +4,67 @@
namespace Dune namespace Dune
{ {
template <class Grid> template <class Grid, class GridReaderImp>
class FileReader class FileReader
{ {
private:
// type of underlying implementation, for internal use only
using Implementation = GridReaderImp;
//! \brief An accessor class to call protected members of reader implementations.
struct Accessor : public Implementation
{
template <class... Args>
static std::unique_ptr<Grid> readImp (Args&&... args)
{
return Implementation::readImpl(std::forward<Args>(args)...);
}
template <class... Args>
static void readFactoryImpl (Args&&... args)
{
return Implementation::readFactoryImpl(std::forward<Args>(args)...);
}
};
public: public:
/// Virtual destructor //! Reads the grid from a file with filename and returns a unique_ptr to the created grid.
virtual ~FileReader() = default; //! Redirects to concrete implementation of derivated class.
template <class... Args>
static std::unique_ptr<Grid> read (const std::string &filename, Args&&... args)
{
return Accessor::readImpl(filename, std::forward<Args>(args)...);
}
//! Reads the grid from a file with filename into a grid-factory.
//! Redirects to concrete implementation of derivated class.
template <class... Args>
static void read (GridFactory<Grid> &factory, const std::string &filename, Args&&... args)
{
Accessor::readFactoryImpl(factory, filename, std::forward<Args>(args)...);
}
/// Reads the grid from a file with `filename` and returns a unique_ptr to the created grid. protected: // default implementations
virtual std::unique_ptr<Grid> read(std::string const& filename)
// Default implementation, redirects to factory read implementation.
template <class... Args>
static std::unique_ptr<Grid> readImpl (const std::string &filename, Args&&... args)
{ {
GridFactory<Grid> factory; GridFactory<Grid> factory;
read(factory, filename); read(factory, filename, std::forward<Args>(args)...);
return std::unique_ptr<Grid>{ factory.createGrid() }; return std::unique_ptr<Grid>{ factory.createGrid() };
} }
/// Reads the grid from a file with `filename` into a grid-factory. // Default implementation for reading into grid-factory: produces a runtime-error.
virtual void read(GridFactory<Grid>& factory, std::string const& filename) = 0; template <class... Args>
static void readFactoryImpl (GridFactory<Grid> &/*factory*/, const std::string &/*filename*/,
Args&&... /*args*/)
{
DUNE_THROW(NotImplemented,
"GridReader using a factory argument not implemented for concrete reader implementation.");
}
}; };
} // end namespace Dune } // end namespace Dune
...@@ -7,6 +7,7 @@ ...@@ -7,6 +7,7 @@
namespace Dune namespace Dune
{ {
/// An abstract base class for LocalFunctions
template <class GridView> template <class GridView>
class VTKLocalFunctionInterface class VTKLocalFunctionInterface
{ {
...@@ -14,16 +15,16 @@ namespace Dune ...@@ -14,16 +15,16 @@ namespace Dune
using Entity = typename GridView::template Codim<0>::Entity; using Entity = typename GridView::template Codim<0>::Entity;
using LocalCoordinate = typename Entity::Geometry::LocalCoordinate; using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
//! Bind the function to the grid entity /// Bind the function to the grid entity
virtual void bind (Entity const& entity) = 0; virtual void bind (Entity const& entity) = 0;
//! Unbind from the currently bound entity /// Unbind from the currently bound entity
virtual void unbind () = 0; virtual void unbind () = 0;
//! Evaluate single component comp in the entity at local coordinates xi /// Evaluate single component comp in the entity at local coordinates xi
virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0; virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0;
//! Virtual destructor /// Virtual destructor
virtual ~VTKLocalFunctionInterface () = default; virtual ~VTKLocalFunctionInterface () = default;
}; };
...@@ -64,14 +65,17 @@ namespace Dune ...@@ -64,14 +65,17 @@ namespace Dune
} }
private: private:
// Evaluate a component of a vector valued data
double evaluateImpl (int comp, LocalCoordinate const& xi, std::true_type) const double evaluateImpl (int comp, LocalCoordinate const& xi, std::true_type) const
{ {
auto y = this->get()(xi); auto y = this->get()(xi);
return y[comp]; return y[comp];
} }
double evaluateImpl (int /*comp*/, LocalCoordinate const& xi, std::false_type) const // Return the scalar values
double evaluateImpl (int comp, LocalCoordinate const& xi, std::false_type) const
{ {
assert(comp == 0);
return this->get()(xi); return this->get()(xi);
} }
}; };
...@@ -97,18 +101,19 @@ namespace Dune ...@@ -97,18 +101,19 @@ namespace Dune
VTKLocalFunction () = default; VTKLocalFunction () = default;
//! Bind the function to the grid entity /// Bind the function to the grid entity
void bind (Entity const& entity) void bind (Entity const& entity)
{ {
this->asInterface().bind(entity); this->asInterface().bind(entity);
} }
//! Unbind from the currently bound entity /// Unbind from the currently bound entity
void unbind () void unbind ()
{ {
this->asInterface().unbind(); this->asInterface().unbind();
} }
/// Evaluate the `comp` component of the Range value at local coordinate `xi`
double evaluate (int comp, LocalCoordinate const& xi) const double evaluate (int comp, LocalCoordinate const& xi) const
{ {
return this->asInterface().evaluate(comp, xi); return this->asInterface().evaluate(comp, xi);
...@@ -117,14 +122,15 @@ namespace Dune ...@@ -117,14 +122,15 @@ namespace Dune
// --------------------------------------------------------------------------- // ---------------------------------------------------------------------------
/// An abstract base class for GlobalFunctions
template <class GridView> template <class GridView>
class VTKFunctionInterface class VTKFunctionInterface
{ {
public: public:
//! Create a local function /// Create a local function
virtual VTKLocalFunction<GridView> makelocalFunction () const = 0; virtual VTKLocalFunction<GridView> makelocalFunction () const = 0;
//! Virtual destructor /// Virtual destructor
virtual ~VTKFunctionInterface () = default; virtual ~VTKFunctionInterface () = default;
}; };
...@@ -156,27 +162,35 @@ namespace Dune ...@@ -156,27 +162,35 @@ namespace Dune
using Super = Functions::TypeErasureBase<VTKFunctionInterface<GridView>, using Super = Functions::TypeErasureBase<VTKFunctionInterface<GridView>,
VTKFunctionImpl<GridView>::template Model>; VTKFunctionImpl<GridView>::template Model>;
template <class F>
using IsGridFunction = decltype(localFunction(std::declval<F>()));
public: public:
template <class F, disableCopyMove<VTKFunction, F> = 0> template <class F, disableCopyMove<VTKFunction, F> = 0>
VTKFunction (F&& f, std::string name, int ncomps = 1) VTKFunction (F&& f, std::string name, int ncomps = 1)
: Super(std::forward<F>(f)) : Super(std::forward<F>(f))
, name_(std::move(name)) , name_(std::move(name))
, ncomps_(ncomps) , ncomps_(ncomps)
{} {
static_assert(Std::is_detected<IsGridFunction,F>::value,
"Requires A GridFunction to be passed to the VTKFunction.");
}
VTKFunction () = default; VTKFunction () = default;
//! Bind the function to the grid entity /// Create a LocalFunction
friend VTKLocalFunction<GridView> localFunction (VTKFunction const& self) friend VTKLocalFunction<GridView> localFunction (VTKFunction const& self)
{ {
return self.asInterface().makelocalFunction(); return self.asInterface().makelocalFunction();
} }
/// Return the number of components of the Range
int ncomps () const int ncomps () const
{ {
return ncomps_; return ncomps_; // TODO: detect automatically. Is this always possible?
} }
/// Return a name associated with the function
std::string name () const std::string name () const
{ {
return name_; return name_;
......
...@@ -9,10 +9,17 @@ ...@@ -9,10 +9,17 @@
namespace Dune namespace Dune
{ {
/// File-Reader for Vtk .vtu files /// File-Reader for Vtk .vtu files
/**
* Reads .vtu files and constructs a grid from the cells stored in the file
* Additionally, stored data can be read.
*
* Assumption on the file structure: Each XML tag must be on a separate line.
**/
template <class Grid> template <class Grid>
class VtkReader class VtkReader
: public FileReader<Grid> : public FileReader<Grid, VtkReader<Grid>>
{ {
// Sections visited during the xml parsing
enum Sections { enum Sections {
NO_SECTION = 0, VTK_FILE, UNSTRUCTURED_GRID, PIECE, POINT_DATA, PD_DATA_ARRAY, CELL_DATA, CD_DATA_ARRAY, NO_SECTION = 0, VTK_FILE, UNSTRUCTURED_GRID, PIECE, POINT_DATA, PD_DATA_ARRAY, CELL_DATA, CD_DATA_ARRAY,
POINTS, POINTS_DATA_ARRAY, CELLS, CELLS_DATA_ARRAY, APPENDED_DATA, XML_NAME, XML_NAME_ASSIGN, XML_VALUE POINTS, POINTS_DATA_ARRAY, CELLS, CELLS_DATA_ARRAY, APPENDED_DATA, XML_NAME, XML_NAME_ASSIGN, XML_VALUE
...@@ -22,42 +29,51 @@ namespace Dune ...@@ -22,42 +29,51 @@ namespace Dune
using GlobalCoordinate = typename Entity::Geometry::GlobalCoordinate; using GlobalCoordinate = typename Entity::Geometry::GlobalCoordinate;
public: public:
/// Constructor /// Constructor. Stores a pointer to the GridFactory.
VtkReader () VtkReader (GridFactory<Grid>& factory)
: buffer_(block_size) : factory_(&factory)
{} {}
/// read the file /// Read the grid from file with `filename` into the GridFactory `factory`
virtual void read (GridFactory<Grid>& factory, std::string const& filename) override; void readFromFile (std::string const& filename);
using FileReader<Grid>::read;
/// Implementation of \ref FileReader interface
static void readFactoryImpl (GridFactory<Grid>& factory, std::string const& filename)
{
VtkReader reader{factory};
reader.readFromFile(filename);
}
private: private:
Sections readCellData (std::ifstream&, // Read values stored on the cells with name `name`
GridFactory<Grid>& /*factory*/, template <class T>
Sections readCellData (std::ifstream& /*input*/,
std::vector<T>& /*values*/,
std::string /*name*/, std::string /*name*/,
Vtk::DataTypes /*type*/, Vtk::DataTypes /*type*/,
std::size_t /*nComponents*/, std::size_t /*nComponents*/,
std::string /*format*/, std::string /*format*/,
std::size_t /*offset*/) std::uint64_t /*offset*/)
{ {
/* does not read anything */ /* does not read anything */
return CD_DATA_ARRAY; return CD_DATA_ARRAY;
} }
Sections readPointData (std::ifstream&, template <class T>
GridFactory<Grid>& /*factory*/, Sections readPointData (std::ifstream& /*input*/,
std::vector<T>& /*values*/,
std::string /*name*/, std::string /*name*/,
Vtk::DataTypes /*type*/, Vtk::DataTypes /*type*/,
std::size_t /*nComponents*/, std::size_t /*nComponents*/,
std::string /*format*/, std::string /*format*/,
std::size_t /*offset*/) std::uint64_t /*offset*/)
{ {
/* does not read anything */ /* does not read anything */
return PD_DATA_ARRAY; return PD_DATA_ARRAY;
} }
// Read vertex coordinates from `input` stream and store in into `factory`
Sections readPoints (std::ifstream& input, Sections readPoints (std::ifstream& input,
GridFactory<Grid>& factory,
std::string name, std::string name,
Vtk::DataTypes type, Vtk::DataTypes type,
std::size_t nComponents, std::size_t nComponents,
...@@ -65,11 +81,10 @@ namespace Dune ...@@ -65,11 +81,10 @@ namespace Dune
std::uint64_t offset); std::uint64_t offset);
template <class T> template <class T>
void readPointsAppended (std::ifstream& input, void readPointsAppended (std::ifstream& input);
GridFactory<Grid>& factory);
// Read cell type, cell offsets and connectivity from `input` stream
Sections readCells (std::ifstream& input, Sections readCells (std::ifstream& input,
GridFactory<Grid>& factory,
std::string name, std::string name,
Vtk::DataTypes type, Vtk::DataTypes type,
std::string format, std::string format,
...@@ -77,13 +92,15 @@ namespace Dune ...@@ -77,13 +92,15 @@ namespace Dune
void readCellsAppended (std::ifstream& input); void readCellsAppended (std::ifstream& input);
// Read data from appended section in vtk file, starting from `offset`
template <class T> template <class T>
void readAppended (std::ifstream& input, std::vector<T>& values, std::uint64_t offset); void readAppended (std::ifstream& input, std::vector<T>& values, std::uint64_t offset);
inline bool isSection (std::string line, // Test whether line belongs to section
std::string key, bool isSection (std::string line,
Sections current, std::string key,
Sections parent = NO_SECTION) Sections current,
Sections parent = NO_SECTION)
{ {
bool result = line.substr(1, key.length()) == key; bool result = line.substr(1, key.length()) == key;
if (result && current != parent) if (result && current != parent)
...@@ -91,27 +108,33 @@ namespace Dune ...@@ -91,27 +108,33 @@ namespace Dune
return result; return result;
} }
std::map<std::string, std::string> parseXml(std::string const& line); // Read attributes from current xml tag
std::map<std::string, std::string> parseXml(std::string const& line, bool& closed);
void createGrid(GridFactory<Grid>& factory) const; // Construct a grid using the GridFactory `factory` and the read vectors
// \ref vec_types, \ref vec_offsets, and \ref vec_connectivity
void createGrid() const;
private: private:
GridFactory<Grid>* factory_;
/// Data format, i.e. ASCII, BINARY or COMPRESSED. Read from xml attributes.
Vtk::FormatTypes format_; Vtk::FormatTypes format_;
std::vector<std::uint8_t> vec_types; // Temporary data to construct the grid elements
std::vector<std::int64_t> vec_offsets; std::vector<std::uint8_t> vec_types; //< VTK cell type ID
std::vector<std::int64_t> vec_connectivity; std::vector<std::int64_t> vec_offsets; //< offset of vertices of cell
std::vector<std::int64_t> vec_connectivity; //< vertex indices of cell
std::size_t numCells_; std::size_t numCells_; //< Number of cells in the grid
std::size_t numVertices_; std::size_t numVertices_; // Number of vertices in the grid
std::size_t numData_;
// offset information for appended data
// map Name -> {DataType,Offset} // map Name -> {DataType,Offset}
std::map<std::string, std::pair<Vtk::DataTypes,std::uint64_t>> offsets_; std::map<std::string, std::pair<Vtk::DataTypes,std::uint64_t>> offsets_;
std::uint64_t offset0_;
std::size_t const block_size = 1024*32; /// Offset of beginning of appended data
std::vector<unsigned char> buffer_; std::uint64_t offset0_;
}; };
} // end namespace Dune } // end namespace Dune
......
...@@ -2,19 +2,20 @@ ...@@ -2,19 +2,20 @@
#include <fstream> #include <fstream>
#include <iterator> #include <iterator>
#include <string> #include <string>
//#include <regex>
#ifdef HAVE_ZLIB #ifdef HAVE_ZLIB
#include <zlib.h> #include <zlib.h>
#endif #endif
#include <dune/common/classname.hh>
#include "utility/filesystem.hh" #include "utility/filesystem.hh"
#include "utility/string.hh" #include "utility/string.hh"
namespace Dune { namespace Dune {
template <class Grid> template <class Grid>
void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filename) void VtkReader<Grid>::readFromFile (std::string const& filename)
{ {
// check whether file exists! // check whether file exists!
if (!filesystem::exists(filename)) if (!filesystem::exists(filename))
...@@ -26,7 +27,7 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -26,7 +27,7 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
data_format = ""; data_format = "";
Vtk::DataTypes data_type = Vtk::UNKNOWN; Vtk::DataTypes data_type = Vtk::UNKNOWN;
std::size_t data_components = 0; std::size_t data_components = 0;
std::size_t data_offset = 0; std::uint64_t data_offset = 0;
std::string file_type = "", std::string file_type = "",
byte_order = "", byte_order = "",
...@@ -40,15 +41,17 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -40,15 +41,17 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
ltrim(line); ltrim(line);
if (isSection(line, "VTKFile", section)) { if (isSection(line, "VTKFile", section)) {
bool closed = false;
auto attr = parseXml(line); auto attr = parseXml(line, closed);
file_type = attr["type"]; file_type = attr["type"];
if (file_type != "UnstructuredGrid") if (file_type != "UnstructuredGrid")
DUNE_THROW(NotImplemented, "Only UnstructuredGrid format implemented. Found: " << file_type); DUNE_THROW(NotImplemented, "Only UnstructuredGrid format implemented. Found: " << file_type);
if (!attr["version"].empty()) if (!attr["version"].empty())
version = std::stod(attr["version"]); assert(std::stod(attr["version"]) == 1.0);
byte_order = attr["byte_order"]; if (!attr["header_type"].empty())
header_type = attr["header_type"]; assert(attr["byte_order"] == "LittleEndian");
if (!attr["header_type"].empty())
assert(attr["header_type"] == "UInt64");
if (!attr["compressor"].empty()) if (!attr["compressor"].empty())
compressor = attr["compressor"]; compressor = attr["compressor"];
section = VTK_FILE; section = VTK_FILE;
...@@ -60,7 +63,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -60,7 +63,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
else if (isSection(line, "/UnstructuredGrid", section, UNSTRUCTURED_GRID)) else if (isSection(line, "/UnstructuredGrid", section, UNSTRUCTURED_GRID))
section = VTK_FILE; section = VTK_FILE;
else if (isSection(line, "Piece", section, UNSTRUCTURED_GRID)) { else if (isSection(line, "Piece", section, UNSTRUCTURED_GRID)) {
auto attr = parseXml(line); bool closed = false;
auto attr = parseXml(line, closed);
numVertices_ = std::stol(attr["NumberOfPoints"]); numVertices_ = std::stol(attr["NumberOfPoints"]);
numCells_ = std::stol(attr["NumberOfCells"]); numCells_ = std::stol(attr["NumberOfCells"]);
section = PIECE; section = PIECE;
...@@ -84,7 +88,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -84,7 +88,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
else if (isSection(line, "/Cells", section, CELLS)) else if (isSection(line, "/Cells", section, CELLS))
section = PIECE; section = PIECE;
else if (line.substr(1,9) == "DataArray") { else if (line.substr(1,9) == "DataArray") {
auto attr = parseXml(line); bool closed = false;
auto attr = parseXml(line, closed);
data_type = Vtk::Map::datatype[attr["type"]]; data_type = Vtk::Map::datatype[attr["type"]];
data_name = attr["Name"]; data_name = attr["Name"];
...@@ -107,7 +112,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -107,7 +112,8 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
assert(data_format == "appended"); assert(data_format == "appended");
} }
if (section == POINT_DATA) if (closed) continue;
else if (section == POINT_DATA)
section = PD_DATA_ARRAY; section = PD_DATA_ARRAY;
else if (section == POINTS) else if (section == POINTS)
section = POINTS_DATA_ARRAY; section = POINTS_DATA_ARRAY;
...@@ -131,37 +137,52 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena ...@@ -131,37 +137,52 @@ void VtkReader<Grid>::read(GridFactory<Grid>& factory, std::string const& filena
DUNE_THROW(Exception, "Wrong section for </DataArray>"); DUNE_THROW(Exception, "Wrong section for </DataArray>");
} }
else if (isSection(line, "AppendedData", section, VTK_FILE)) { else if (isSection(line, "AppendedData", section, VTK_FILE)) {
auto attr = parseXml(line); bool closed = false;