FactoryParametrization.hpp 2.82 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
#pragma once

#include <type_traits>
#include <vector>

#include <dune/grid/common/gridfactory.hh>
#include <amdis/Output.hpp>

namespace AMDiS
{
  template <class GridType, class Projection>
  struct Parametrization
  {
    using Grid = GridType;
    using ctype = typename Grid::ctype;

    static constexpr int dimension = Grid::dimension;
    static constexpr int dimensionworld = Grid::dimensionworld;

    using LeafIntersection = typename Grid::LeafIntersection;

    template <int cd>
    using Codim = typename Grid::template Codim<cd>;
  };

} // end namespace AMDiS


namespace Dune
{
  /// \brief A factory class for \ref SimpleGrid.
  template <class Grid, class Projection>
  class GridFactory<AMDiS::Parametrization<Grid, Projection>>
      : public GridFactoryInterface<AMDiS::Parametrization<Grid, Projection>>
  {
    using GridWrapper = AMDiS::Parametrization<Grid, Projection>;

    /// Type used by the grid for coordinates
    using ctype = typename Grid::ctype;

    static constexpr int dim = Grid::dimension;
    static constexpr int dow = Grid::dimensionworld;

    using LocalCoordinate = FieldVector<ctype, dim>;
    using GlobalCoordinate = FieldVector<ctype, dow>;

    using ProjectionBase = VirtualFunction<LocalCoordinate, GlobalCoordinate>;

  public:

    /// Default constructor
52
    GridFactory() = default;
Praetorius, Simon's avatar
Praetorius, Simon committed
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73

    /// Insert a vertex into the coarse grid
    virtual void insertVertex(GlobalCoordinate const& pos) override
    {
      GlobalCoordinate new_pos(pos);
      Projection::project(new_pos);

      coordinates.push_back(new_pos);
      factory.insertVertex(new_pos);
    }

    /// Insert an element into the coarse grid.
    virtual void insertElement(GeometryType const& type,
                               std::vector<unsigned int> const& vertices) override
    {
      assert( type.isSimplex() );

      std::vector<GlobalCoordinate> corners(vertices.size());
      for (std::size_t i = 0; i < vertices.size(); ++i)
        corners[i] = coordinates[vertices[i]];

74
      factory.insertElement(type, vertices, std::make_shared<Projection>(std::move(corners)) );
Praetorius, Simon's avatar
Praetorius, Simon committed
75
76
77
78
79
80
81
82
    }

    virtual void insertBoundarySegment(std::vector<unsigned int> const& /*vertices*/) override
    {
      AMDiS::warning("Not implemented!");
    }

    /// Finalize grid creation and hand over the grid.
83
    std::unique_ptr<Grid> create()
Praetorius, Simon's avatar
Praetorius, Simon committed
84
    {
85
      return std::unique_ptr<Grid>(factory.createGrid());
Praetorius, Simon's avatar
Praetorius, Simon committed
86
87
    }

88
89
90
#if DUNE_VERSION_GT(DUNE_GRID,2,6)
    virtual ToUniquePtr<GridWrapper> createGrid() override
#else
Praetorius, Simon's avatar
Praetorius, Simon committed
91
    virtual GridWrapper* createGrid() override
92
#endif
Praetorius, Simon's avatar
Praetorius, Simon committed
93
94
95
96
97
98
99
100
101
102
103
104
    {
      AMDiS::warning("Should not be created. Use non-virtual method `create()` instead, to create the underlying grid!");
      return nullptr;
    }

  private:
    // buffers for the mesh data
    std::vector<GlobalCoordinate> coordinates;
    GridFactory<Grid> factory;
  };

} // end namespace Dune