DefaultGridView.hpp 9.11 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
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
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
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
#pragma once

#include <dune/common/typetraits.hh>
#include <dune/common/exceptions.hh>

#include <dune/grid/common/capabilities.hh>
#include <dune/grid/common/gridview.hh>

namespace AMDiS
{
  /*
   * NOTE: this is a modification of dune/grid/common/defaultgridview.hh since we
   * want all implementation specific methods to be put into the grid class and not
   * into the entity class. See ibegin(), iend().
   */

  template <class GridImp>
  class DefaultLevelGridView;

  template <class GridImp>
  class DefaultLeafGridView;


  template <class GridImp>
  struct DefaultLevelGridViewTraits
  {
    using GridViewImp = DefaultLevelGridView<GridImp>;

    /// type of the grid
    using Grid = std::remove_const_t<GridImp>;

    /// type of the index set
    using IndexSet = typename Grid::Traits::LevelIndexSet;

    /// type of the intersection
    using Intersection = typename Grid::Traits::LevelIntersection;

    /// type of the intersection iterator
    using IntersectionIterator = typename Grid::Traits::LevelIntersectionIterator;

    /// type of the collective communication
    using CollectiveCommunication = typename Grid::Traits::CollectiveCommunication;

    template <int cd>
    struct Codim
    {
      using Entity = typename Grid::Traits::template Codim<cd>::Entity;
      using Geometry = typename Grid::template Codim<cd>::Geometry;
      using LocalGeometry = typename Grid::template Codim<cd>::LocalGeometry;

      /// Define types needed to iterate over entities of a given partition type
      template <Dune::PartitionIteratorType pit>
      struct Partition
      {
        /// iterator over a given codim and partition type
        using Iterator = typename Grid::template Codim<cd>::template Partition<pit>::LevelIterator;
      };

      using Iterator = typename Partition<Dune::All_Partition>::Iterator;
    };

    enum { conforming = Dune::Capabilities::isLevelwiseConforming<Grid>::v };
  };


  template <class GridImp>
  class DefaultLevelGridView
  {
  public:
    using Traits = DefaultLevelGridViewTraits<GridImp>;
    using Grid = typename Traits::Grid;
    using IndexSet = typename Traits::IndexSet;
    using Intersection = typename Traits::Intersection;
    using IntersectionIterator = typename Traits::IntersectionIterator;
    using CollectiveCommunication = typename Traits::CollectiveCommunication;

    template <int cd>
    using Codim = typename Traits::template Codim<cd>;

    enum { conforming = Traits::conforming };

  public:
    DefaultLevelGridView(Grid const& grid, int level)
      : grid_(&grid)
      , level_(level)
    {}

    /// Obtain a const reference to the underlying hierarchic grid
    Grid const& grid() const
    {
      assert(grid_);
      return *grid_;
    }

    /// Obtain the index set
    IndexSet const& indexSet() const
    {
      return grid().levelIndexSet(level_);
    }

    /// Obtain number of entities in a given codimension
    int size(int codim) const
    {
      return grid().size(level_, codim);
    }

    /// Obtain number of entities with a given geometry type
    int size(Dune::GeometryType const& type) const
    {
      return grid().size(level_, type);
    }

    /// Obtain begin iterator for this view
    template <int cd, Dune::PartitionIteratorType pit = Dune::All_Partition>
    typename Codim<cd>::template Partition<pit>::Iterator begin() const
    {
      return grid().template lbegin<cd, pit>(level_);
    }

    /// Obtain end iterator for this view
    template <int cd, Dune::PartitionIteratorType pit = Dune::All_Partition>
    typename Codim<cd>::template Partition<pit>::Iterator end() const
    {
      return grid().template lend<cd, pit>(level_);
    }

    /// Obtain begin intersection iterator with respect to this view
    IntersectionIterator ibegin(typename Codim<0>::Entity const& entity) const
    {
      return grid().ilevelbegin(entity);
    }

    /// Obtain end intersection iterator with respect to this view
    IntersectionIterator iend(typename Codim<0>::Entity const& entity) const
    {
      return grid().ilevelend(entity);
    }

    /// Obtain collective communication object
    CollectiveCommunication const& comm() const
    {
      return grid().comm();
    }

    /// Return size of the overlap region for a given codim on the grid view.
    int overlapSize(int codim) const
    {
      return grid().overlapSize(level_, codim);
    }

    /// Return size of the ghost region for a given codim on the grid view.
    int ghostSize(int codim) const
    {
      return grid().ghostSize(level_, codim);
    }

    /// Communicate data on this view
    template <class DataHandleImp, class DataType>
    void communicate(Dune::CommDataHandleIF<DataHandleImp, DataType>& data,
                     Dune::InterfaceType iftype,
                     Dune::CommunicationDirection dir) const
    {
      return grid().communicate(data, iftype, dir, level_);
    }

  private:
    Grid const* grid_;
    int level_;
  };


  template <class GridImp>
  struct DefaultLeafGridViewTraits
  {
    using GridViewImp = DefaultLeafGridView<GridImp>;

    /// type of the grid
    using Grid = std::remove_const_t<GridImp>;

    /// type of the index set
    using IndexSet = typename Grid::Traits::LeafIndexSet;

    /// type of the intersection
    using Intersection = typename Grid::Traits::LeafIntersection;

    /// type of the intersection iterator
    using IntersectionIterator = typename Grid::Traits::LeafIntersectionIterator;

    /// type of the collective communication
    using CollectiveCommunication = typename Grid::Traits::CollectiveCommunication;

    template <int cd>
    struct Codim
    {
      using Entity = typename Grid::Traits::template Codim<cd>::Entity;
      using Geometry = typename Grid::template Codim<cd>::Geometry;
      using LocalGeometry = typename Grid::template Codim<cd>::LocalGeometry;

      /// Define types needed to iterate over entities of a given partition type
      template <Dune::PartitionIteratorType pit>
      struct Partition
      {
        /// iterator over a given codim and partition type
        using Iterator = typename Grid::template Codim<cd>::template Partition<pit>::LeafIterator;
      };

      using Iterator = typename Partition<Dune::All_Partition>::Iterator;
    };

    enum { conforming = Dune::Capabilities::isLeafwiseConforming<Grid>::v };
  };


  template <class GridImp>
  class DefaultLeafGridView
  {
  public:
    using Traits = DefaultLeafGridViewTraits<GridImp>;
    using Grid = typename Traits::Grid;
    using IndexSet = typename Traits::IndexSet;
    using Intersection = typename Traits::Intersection;
    using IntersectionIterator = typename Traits::IntersectionIterator;
    using CollectiveCommunication = typename Traits::CollectiveCommunication;

    template <int cd>
    using Codim = typename Traits::template Codim<cd>;

    enum { conforming = Traits::conforming };

  public:
    DefaultLeafGridView(Grid const& grid)
      : grid_(&grid)
    {}

    /// obtain a const reference to the underlying hierarchic grid
    Grid const& grid() const
    {
      assert(grid_);
      return *grid_;
    }

    /// Obtain the index set
    IndexSet const& indexSet() const
    {
      return grid().leafIndexSet();
    }

    /// Obtain number of entities in a given codimension
    int size(int codim) const
    {
      return grid().size(codim);
    }

    /// Obtain number of entities with a given geometry type
    int size(Dune::GeometryType const& type) const
    {
      return grid().size(type);
    }


    /// Obtain begin iterator for this view
    template <int cd, Dune::PartitionIteratorType pit = Dune::All_Partition>
    typename Codim<cd>::template Partition<pit>::Iterator begin() const
    {
      return grid().template leafbegin<cd, pit>();
    }

    /// Obtain end iterator for this view
    template <int cd, Dune::PartitionIteratorType pit = Dune::All_Partition>
    typename Codim<cd>::template Partition<pit>::Iterator end() const
    {
      return grid().template leafend<cd, pit>();
    }

    /// Obtain begin intersection iterator with respect to this view
    IntersectionIterator ibegin(typename Codim<0>::Entity const& entity) const
    {
      return grid().ileafbegin(entity);
    }

    /// Obtain end intersection iterator with respect to this view
    IntersectionIterator iend(typename Codim<0>::Entity const& entity) const
    {
      return grid().ileafend(entity);
    }

    /// Obtain collective communication object
    CollectiveCommunication const& comm() const
    {
      return grid().comm();
    }

    /// Return size of the overlap region for a given codim on the grid view.
    int overlapSize(int codim) const
    {
      return grid().overlapSize(codim);
    }

    /// Return size of the ghost region for a given codim on the grid view.
    int ghostSize(int codim) const
    {
      return grid().ghostSize(codim);
    }

    /// Communicate data on this view
    template <class DataHandleImp, class DataType>
    void communicate(Dune::CommDataHandleIF<DataHandleImp, DataType>& data,
                     Dune::InterfaceType iftype,
                     Dune::CommunicationDirection dir) const
    {
      return grid().communicate(data, iftype, dir);
    }

  private:
    Grid const* grid_;
  };

} // end namespace AMDiS