MacroGridFactory.hpp 5.54 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
#pragma once

#include <dune/common/typeutilities.hh>
#include <dune/common/version.hh>
#include <dune/grid/utility/structuredgridfactory.hh>

#if ! DUNE_VERSION_GT(DUNE_GRID,2,6)
  #include <dune/grid/yaspgrid.hh>
#endif

namespace Dune
{
  template <class GridType>
  struct CubedWrapper : public GridType {};


  template <class GridType>
  class GridFactory<CubedWrapper<GridType>>
      : public GridFactoryInterface<GridType>
  {
    using Self = GridFactory;
22
    using Super = GridFactoryInterface<GridType>;
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
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
    using ctype = typename GridType::ctype;

    enum { dim = GridType::dimension };
    enum { dimworld = GridType::dimensionworld };

  public:
    template <class... Args, disableCopyMove<Self, Args...> = 0>
    GridFactory (Args&&... args)
      : factory_(std::make_shared<GridFactory<GridType>>(std::forward<Args>(args)...))
    {}

    GridFactory (GridFactory<GridType>& factory)
      : factory_(stackobject_to_shared_ptr(factory))
    {}

    /// \brief Insert a vertex into the coarse grid
    void insertVertex (const FieldVector<ctype,dimworld>& pos) override
    {
      factory_->insertVertex(pos);
    }

    /// \brief Insert simplex elements into the coarse grid
    /**
     * Creates a simplex subdividion of the cube element.
     *
     * \param type      The GeometryType of the box grid
     * \param vertices  Indices of the cube corners
     **/
    void insertElement (const GeometryType& type,
                        const std::vector<unsigned int>& vertices) override
    {
      // triangulation of reference cube
      static const auto reference_cubes = std::make_tuple(
        std::array<std::array<int,2>, 1>{std::array<int,2>{0,1}},
        std::array<std::array<int,3>, 2>{std::array<int,3>{3,0,1}, std::array<int,3>{0,3,2}},
        std::array<std::array<int,4>, 6>{std::array<int,4>{0,7,3,1}, std::array<int,4>{0,7,5,1},
                                         std::array<int,4>{0,7,5,4}, std::array<int,4>{0,7,3,2},
                                         std::array<int,4>{0,7,6,2}, std::array<int,4>{0,7,6,4}} );

      assert(type == GeometryTypes::cube(dim));

      auto const& simplices = std::get<dim-1>(reference_cubes);
      thread_local std::vector<unsigned int> corners(dim+1);
      for (auto const& simplex : simplices) {
        for (std::size_t i = 0; i < simplex.size(); ++i)
          corners[i] = vertices[simplex[i]];

        factory_->insertElement(GeometryTypes::simplex(dim), corners);
      }
    }

74
75
    using Super::insertElement;

76
77
78
79
80
81
82
    /// \brief insert a boundary segment
    // TODO: maybe split boundary segment in simplices
    void insertBoundarySegment (const std::vector<unsigned int>& vertices) override
    {
      factory_->insertBoundarySegment(vertices);
    }

83
84
    using Super::insertBoundarySegment;

85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
    /// \brief Finalize grid creation and hand over the grid
#if DUNE_VERSION_GT(DUNE_GRID,2,6)
    ToUniquePtr<GridType> createGrid () override
#else
    GridType* createGrid () override
#endif
    {
      return factory_->createGrid();
    }

  private:
    std::shared_ptr<GridFactory<GridType>> factory_;
  };

} // end namespace Dune


namespace AMDiS
{
  template <class GridType>
  class MacroGridFactory
  {
    using ctype = typename GridType::ctype;

    enum { dim = GridType::dimension };
    enum { dimworld = GridType::dimensionworld };

  public:
#if DUNE_VERSION_GT(DUNE_GRID,2,6)
    /// \brief insert structured simplex grid into grid factory
    static std::unique_ptr<GridType> createSimplexGrid (Dune::GridFactory<GridType>& originalFactory,
                                                        const Dune::FieldVector<ctype,dimworld>& lowerLeft,
                                                        const Dune::FieldVector<ctype,dimworld>& upperRight,
                                                        const std::array<unsigned int,dim>& numElements)
    {
      Dune::GridFactory<Dune::CubedWrapper<GridType>> factory(originalFactory);
      Dune::StructuredGridFactory<Dune::CubedWrapper<GridType>>::createCubeGrid(factory, lowerLeft, upperRight, numElements);
      return std::unique_ptr<GridType>(factory.createGrid());
    }
#endif

    /// \brief Create a structured simplex grid
    static std::unique_ptr<GridType> createSimplexGrid (const Dune::FieldVector<ctype,dimworld>& lowerLeft,
                                                        const Dune::FieldVector<ctype,dimworld>& upperRight,
                                                        const std::array<unsigned int,dim>& numElements)
    {
      Dune::GridFactory<Dune::CubedWrapper<GridType>> factory;
#if DUNE_VERSION_GT(DUNE_GRID,2,6)
      Dune::StructuredGridFactory<Dune::CubedWrapper<GridType>>::createCubeGrid(factory, lowerLeft, upperRight, numElements);
#else
      // fallback implementation using temporary YaspGrid
      using TempGrid = Dune::YaspGrid<dim,Dune::EquidistantOffsetCoordinates<double,dim>>;
      auto grid = Dune::StructuredGridFactory<TempGrid>::createCubeGrid(lowerLeft, upperRight, numElements);
      for (auto const& v : vertices(grid->leafGridView()))
        factory.insertVertex(v.geometry().corner(0));

      auto const& indexSet = grid->leafIndexSet();
      for (auto const& e : elements(grid->leafGridView()))
      {
        thread_local std::vector<unsigned int> vertices;
        vertices.resize(e.subEntities(dim));
        for (unsigned int i = 0; i < e.subEntities(dim); ++i)
          vertices[i] = indexSet.subIndex(e,i,dim);
        factory.insertElement(Dune::GeometryTypes::cube(dim), vertices);
      }
#endif
      return std::unique_ptr<GridType>(factory.createGrid());
    }
  };

} // end namespace AMDiS