dune-vtk 2.8
Loading...
Searching...
No Matches
continuousgridcreator.hh
Go to the documentation of this file.
1#pragma once
2
3#include <cassert>
4#include <cstdint>
5#include <limits>
6#include <vector>
7
8#include <dune/common/exceptions.hh>
9#include <dune/common/hybridutilities.hh>
10#include <dune/grid/common/gridfactory.hh>
11
12#include <dune/vtk/types.hh>
15
16namespace Dune
17{
18 namespace Vtk
19 {
20 // Create a grid where the input points and connectivity is already
21 // connected correctly.
22 template <class Grid>
24 : public GridCreatorInterface<Grid, ContinuousGridCreator<Grid>>
25 {
29
30 using Super::Super;
31 using Super::factory;
32
33 void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
34 std::vector<std::uint64_t> const& /*point_ids*/)
35 {
36 for (auto const& p : points)
37 factory().insertVertex(p);
38 }
39
40 void insertElementsImpl (std::vector<std::uint8_t> const& types,
41 std::vector<std::int64_t> const& offsets,
42 std::vector<std::int64_t> const& connectivity)
43 {
44 std::size_t idx = 0;
45 for (std::size_t i = 0; i < types.size(); ++i) {
46 auto type = Vtk::to_geometry(types[i]);
47 Vtk::CellType cellType{type};
48 [[maybe_unused]] auto refElem = referenceElement<double,Grid::dimension>(type);
49
50 int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
51 assert(nNodes == refElem.size(Grid::dimension));
52 std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
53 for (int j = 0; j < nNodes; ++j)
54 vtk_cell.push_back( connectivity[idx++] );
55
56 if (cellType.noPermutation())
57 factory().insertElement(type,vtk_cell);
58 else {
59 // apply index permutation
60 std::vector<unsigned int> cell(nNodes);
61 for (int j = 0; j < nNodes; ++j)
62 cell[j] = vtk_cell[cellType.permutation(j)];
63
64 factory().insertElement(type,cell);
65 }
66 }
67 }
68 };
69
70 // deduction guides
71 template <class Grid>
72 ContinuousGridCreator(GridFactory<Grid>&)
74
75 template <class GridType, class FieldType, class Context>
76 struct AssociatedGridFunction<ContinuousGridCreator<GridType>, FieldType, Context>
77 {
79 };
80
81 } // end namespace Vtk
82} // end namespace Dune
Definition: writer.hh:13
GeometryType to_geometry(std::uint8_t cell)
Definition: types.cc:146
Base class for grid creators in a CRTP style.
Definition: gridcreatorinterface.hh:25
GridFactory< Grid > & factory()
Return the associated GridFactory.
Definition: gridcreatorinterface.hh:77
typename Grid::template Codim< 0 >::Entity::Geometry::GlobalCoordinate GlobalCoordinate
Definition: gridcreatorinterface.hh:28
Definition: continuousgridcreator.hh:25
GridFactory< Grid > & factory()
Return the associated GridFactory.
Definition: gridcreatorinterface.hh:77
void insertElementsImpl(std::vector< std::uint8_t > const &types, std::vector< std::int64_t > const &offsets, std::vector< std::int64_t > const &connectivity)
Definition: continuousgridcreator.hh:40
typename Super::GlobalCoordinate GlobalCoordinate
Definition: continuousgridcreator.hh:28
void insertVerticesImpl(std::vector< GlobalCoordinate > const &points, std::vector< std::uint64_t > const &)
Definition: continuousgridcreator.hh:33
Type-Traits to associate a GridFunction to a GridCreator.
Definition: gridfunctions/common.hh:26
A GridFunction representing data stored on the grid vertices in a continuous manner.
Definition: continuousgridfunction.hh:16
Mapping of Dune geometry types to VTK cell types.
Definition: types.hh:160