8 #include <dune/common/shared_ptr.hh>
9 #include <dune/common/typelist.hh>
10 #include <dune/common/typeutilities.hh>
34 template <
class Gr
id,
class GC = Vtk::ContinuousGr
idCreator<Gr
id>,
class FieldType =
double>
40 NO_SECTION = 0, VTK_FILE, UNSTRUCTURED_GRID, PIECE, POINT_DATA, PD_DATA_ARRAY, CELL_DATA, CD_DATA_ARRAY,
41 POINTS, POINTS_DATA_ARRAY, CELLS, CELLS_DATA_ARRAY, APPENDED_DATA, XML_NAME, XML_NAME_ASSIGN, XML_VALUE
45 struct DataArrayAttributes
49 unsigned int components = 1;
50 std::uint64_t offset = 0;
51 Sections section = NO_SECTION;
55 using GlobalCoordinate =
typename GC::GlobalCoordinate;
61 template <
class Context>
81 template <
class... Args,
82 std::enable_if_t<std::is_constructible<
GridCreator, Args...>::value,
int> = 0>
89 :
VtkReader(stackobject_to_shared_ptr(creator))
93 explicit VtkReader (std::shared_ptr<GridCreator> creator)
94 : creator_(std::move(creator))
106 void read (std::string
const& filename,
bool fillCreator =
true);
118 return creator_->createGrid();
123 GridFunction<Vtk::PointContext>
getPointData (std::string
const& name)
const
125 auto data_it = dataArray_.find(
"PointData." + name);
126 auto point_it = pointData_.find(
"PointData." + name);
127 VTK_ASSERT_MSG(data_it != dataArray_.end() && point_it != pointData_.end(),
128 "The data to extract is not found in point-data. Try `getCellData()` instead!");
129 VTK_ASSERT(data_it->second.section == POINT_DATA);
131 return {*creator_, point_it->second, data_it->second.components,
132 vec_types, vec_offsets, vec_connectivity};
139 std::vector<DataArrayAttributes> attributes;
140 attributes.reserve(pointData_.size());
141 for (
auto const& da : dataArray_) {
142 if (da.second.section == POINT_DATA)
143 attributes.push_back(da.second);
150 GridFunction<Vtk::CellContext>
getCellData (std::string
const& name)
const
152 auto data_it = dataArray_.find(
"CellData." + name);
153 auto cell_it = cellData_.find(
"CellData." + name);
154 VTK_ASSERT_MSG(data_it != dataArray_.end() && cell_it != cellData_.end(),
155 "The data to extract is not found in cell-data. Try `getPointData()` instead!");
156 VTK_ASSERT(data_it->second.section == CELL_DATA);
158 return {*creator_, cell_it->second, data_it->second.components,
159 vec_types, vec_offsets, vec_connectivity};
166 std::vector<DataArrayAttributes> attributes;
167 attributes.reserve(cellData_.size());
168 for (
auto const& da : dataArray_) {
169 if (da.second.section == CELL_DATA)
170 attributes.push_back(da.second);
199 std::vector<std::string>
const&
pieces ()
const
206 static void fillFactoryImpl (GridFactory<Grid>& factory, std::string
const& filename)
209 reader.
read(filename);
215 Sections readCellData (std::ifstream& input, std::string
id);
217 template <
class F,
class H>
218 void readCellDataAppended (MetaType<F>, MetaType<H>, std::ifstream& input, std::string
id);
221 Sections readPointData (std::ifstream& input, std::string
id);
223 template <
class F,
class H>
224 void readPointDataAppended (MetaType<F>, MetaType<H>, std::ifstream& input, std::string
id);
228 Sections readPoints (std::ifstream& input, std::string
id);
230 template <
class F,
class H>
231 void readPointsAppended (MetaType<F>, MetaType<H>, std::ifstream& input);
235 Sections readCells (std::ifstream& input, std::string
id);
238 void readCellsAppended (MetaType<H>, std::ifstream& input);
241 template <
class FloatType,
class HeaderType>
242 void readAppended (std::ifstream& input, std::vector<FloatType>& values, HeaderType offset);
245 bool isSection (std::string line,
248 Sections parent = NO_SECTION)
const
250 bool result = line.substr(1, key.length()) == key;
251 if (result && current != parent)
252 DUNE_THROW(Exception ,
"<" << key <<
"> in wrong section." );
257 std::string toString (Sections s)
const;
260 std::uint64_t findAppendedDataPosition (std::ifstream& input)
const;
263 std::map<std::string, std::string> parseXml (std::string
const& line,
bool& closed);
270 return MPIHelper::getCollectiveCommunication();
274 std::shared_ptr<GridCreator> creator_;
283 std::vector<GlobalCoordinate> vec_points;
284 std::vector<std::uint64_t> vec_point_ids;
285 std::vector<std::uint8_t> vec_types;
286 std::vector<std::int64_t> vec_offsets;
287 std::vector<std::int64_t> vec_connectivity;
289 std::size_t numberOfCells_ = 0;
290 std::size_t numberOfPoints_ = 0;
294 std::map<std::string, DataArrayAttributes> dataArray_;
298 std::map<std::string, std::vector<FieldType>> pointData_;
299 std::map<std::string, std::vector<FieldType>> cellData_;
302 std::vector<std::string> pieces_;
305 std::uint64_t offset0_ = 0;
311 template <
class Gr
id>
315 template <
class GridCreator,
316 class = std::void_t<typename GridCreator::Grid>>
320 template <
class GridCreator,
321 class = std::void_t<typename GridCreator::Grid>>
Macro for wrapping error checks and throwing exceptions.
#define VTK_ASSERT_MSG(cond, text)
check if condition cond holds; otherwise, throw a VtkError with a message.
Definition: errors.hh:19
#define VTK_ASSERT(cond)
check if condition cond holds; otherwise, throw a VtkError.
Definition: errors.hh:29
Definition: datacollectorinterface.hh:9
VtkReader(GridFactory< Grid > &) -> VtkReader< Grid, Vtk::ContinuousGridCreator< Grid >>
FormatTypes
Type used for representing the output format.
Definition: types.hh:21
CompressorTypes
Definition: types.hh:139
DataTypes
Definition: types.hh:50
Definition: filereader.hh:16
static void fillFactoryImpl(GridFactory< Grid > &, const std::string &, Args &&...)
Definition: filereader.hh:71
void type
Definition: gridfunctions/common.hh:27
File-Reader for Vtk unstructured .vtu files.
Definition: vtkreader.hh:37
VtkReader(Args &&... args)
Constructor. Creates a new GridCreator with the passed factory.
Definition: vtkreader.hh:83
void read(std::string const &filename, bool fillCreator=true)
Read the grid from file with filename into the GridCreator.
Definition: vtkreader.impl.hh:20
VtkReader(std::shared_ptr< GridCreator > creator)
Constructor. Stores the shared_ptr.
Definition: vtkreader.hh:93
void readParallelFileFromStream(std::ifstream &input, int rank, int size, bool create=true)
Read the grid from and input stream, referring to a .pvtu file, into the GridFactory factory_.
Definition: vtkreader.impl.hh:249
std::unique_ptr< Grid > createGrid() const
Definition: vtkreader.hh:116
VtkReader(GridCreator &creator)
Constructor. Stores the references in a non-destroying shared_ptr.
Definition: vtkreader.hh:88
GC GridCreator
Definition: vtkreader.hh:65
std::vector< DataArrayAttributes > getPointDataAttributes() const
Definition: vtkreader.hh:137
GridFunction< Vtk::CellContext > getCellData(std::string const &name) const
Definition: vtkreader.hh:150
void fillGridCreator(bool insertPieces=true)
Definition: vtkreader.impl.hh:681
GridFunction< Vtk::PointContext > PointGridFunction
GridFunction representing the data stored on the points in the file.
Definition: vtkreader.hh:68
GridFunction< Vtk::PointContext > getPointData(std::string const &name) const
Definition: vtkreader.hh:123
void readSerialFileFromStream(std::ifstream &input, bool create=true)
Read the grid from an input stream, referring to a .vtu file, into the GridFactory factory_.
Definition: vtkreader.impl.hh:43
std::vector< DataArrayAttributes > getCellDataAttributes() const
Definition: vtkreader.hh:164
GridFunction< Vtk::CellContext > CellGridFunction
GridFunction representing the data stored on the cells in the file.
Definition: vtkreader.hh:71
std::vector< std::string > const & pieces() const
Return the filenames of parallel pieces.
Definition: vtkreader.hh:199
GridCreator & gridCreator()
Obtains the creator of the reader.
Definition: vtkreader.hh:109