BoundaryManager.hpp 5.11 KB
Newer Older
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
#pragma once

#include <memory>
#include <vector>

#include <dune/common/hybridutilities.hh>
#include <dune/common/std/type_traits.hh>

#include <amdis/Boundary.hpp>
#include <amdis/common/Concepts.hpp>

namespace AMDiS
{
  class BoundaryManagerBase
  {
  public:
    BoundaryManagerBase(std::size_t numBoundarySegments)
      : boundaryIds_(numBoundarySegments, BoundaryType{1})
    {}

    /// Return the stored boundary id for the given intersection
    template <class Intersection>
    BoundaryType boundaryId(Intersection const& intersection) const
    {
      return boundaryIds_[intersection.boundarySegmentIndex()];
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
28
29
30
31
32
    std::vector<BoundaryType> const boundaryIds() const
    {
      return boundaryIds_;
    }

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
  protected:
    std::vector<BoundaryType> boundaryIds_; // maps a boundarySegementIndex to an ID
  };


  /// Manage boundary ids of boundary segments in a grid
  /**
   * Manager for bounary IDs, that can be initialized from different sources
   * - cube domains can be assigned box boundaries, by specifying left-right,
   *   front-back, and bottom-top ids
   * - An indicator from global coodinates that returns an id. The indicator function
   *   is evaluated in the center points of intersections.
   * - A predicate that return true/false for a boundary part and sets the given id
   * - Read ids from the grid. Those may be initialized by some grid readers. Note:
   *   not all grids support `intersection.boundaryId()`, but e.g. AlbertaGrid and
   *   ALUGrid do support this. Can be read by DGF parser and AlbertaGrid constructor
   *   from macroFileName.
   **/
  template <class G>
  class BoundaryManager
      : public BoundaryManagerBase
  {
    using Super = BoundaryManagerBase;

  public:
    using Grid = G;
    enum { dim = Grid::dimension };
    enum { dow = Grid::dimensionworld };

    using Segment = typename Grid::LevelGridView::Intersection;
    using Domain = typename Grid::template Codim<0>::Geometry::GlobalCoordinate;

  public:
    /// Constructor. Stores a shared pointer to the grid and initializes all
    /// boundary IDs to the value stored in the grid
    BoundaryManager(std::shared_ptr<G> const& grid)
      : Super(grid->numBoundarySegments())
      , grid_(grid)
    {
      setBoundaryId();
    }

    /// Set boundary ids [left,right, front,back, bottom,top] for cube domains
    void setBoxBoundary(std::array<BoundaryType, 2*dow> const& ids)
    {
78
      auto gv = grid_->leafGridView();
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
      for (auto const& e : elements(gv))
      {
        for (auto const& segment : intersections(gv,e)) {
          if (!segment.boundary())
            continue;

          auto n = segment.centerUnitOuterNormal();
          auto index = segment.boundarySegmentIndex();

          for (int i = 0; i < dow; ++i) {
            if (n[i] < -0.5)
              boundaryIds_[index] = ids[2*i];
            else if (n[i] > 0.5)
              boundaryIds_[index] = ids[2*i+1];
          }
        }
      }
    }


    /// Set indicator(center) for all boundary intersections
    template <class Indicator,
      REQUIRES(Concepts::Functor<Indicator, int(Domain)>) >
    void setIndicator(Indicator const& indicator)
    {
104
      auto gv = grid_->leafGridView();
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
      for (auto const& e : elements(gv))
      {
        for (auto const& segment : intersections(gv,e)) {
          if (!segment.boundary())
            continue;

          auto index = segment.boundarySegmentIndex();
          boundaryIds_[index] = indicator(segment.geometry().center());
        }
      }
    }


    /// Set id for all boundary intersections with pred(center) == true
    template <class Predicate,
      REQUIRES(Concepts::Functor<Predicate, bool(Domain)>) >
    void setPredicate(Predicate const& pred, BoundaryType id)
    {
123
      auto gv = grid_->leafGridView();
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
      for (auto const& e : elements(gv))
      {
        for (auto const& segment : intersections(gv,e)) {
          if (!segment.boundary())
            continue;

          auto index = segment.boundarySegmentIndex();
          if (pred(segment.geometry().center()))
            boundaryIds_[index] = id;
        }
      }
    }


    template <class I>
    using HasBoundaryId = decltype(std::declval<I>().boundaryId());

    /// Set boundary ids as stored in the grid, e.g. read by grid reader
    void setBoundaryId()
    {
      if (!Dune::Std::is_detected<HasBoundaryId, Segment>::value)
        return;

147
      auto gv = grid_->leafGridView();
148
149
150
151
152
153
      for (auto const& e : elements(gv))
      {
        for (auto const& segment : intersections(gv,e)) {
          if (!segment.boundary())
            continue;

154
          Dune::Hybrid::ifElse(Dune::Std::is_detected<HasBoundaryId, Segment>{}, [&](auto id) -> void {
155
156
157
            auto index = segment.boundarySegmentIndex();
            boundaryIds_[index] = id(segment).boundaryId();
          });
158
159
160
161
        }
      }
    }

162
163
164
165
166
167
168
    template <class I>
    void setBoundaryIds(std::vector<I> const& ids)
    {
      test_exit(ids.size() == boundaryIds_.size(), "Number of boundary IDs does not match!");
      std::copy(ids.begin(), ids.end(), boundaryIds_.begin());
    }

169
170
171
172
173
174
  private:
    std::shared_ptr<Grid> grid_;
    using Super::boundaryIds_;
  };

} // end namespace AMDiS