Marker.hpp 8.06 KB
Newer Older
1
#pragma once
2

3
4
5
6
7
8
9
10
11
12
#include <string>

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

#include <amdis/AdaptInfo.hpp>
#include <amdis/Flag.hpp>
#include <amdis/Initfile.hpp>

namespace AMDiS {

13

14
15
16
17
18
19
  /**
   * \ingroup Adaption
   *
   * \brief
   * Base class for all scalar markers.
   */
20
  template <class Traits>
21
22
  class Marker
  {
23
  protected:
24
25
26
    using GlobalBasis = typename Traits::GlobalBasis;
    using GridView = typename GlobalBasis::GridView;
    using Grid     = typename GridView::Grid;
27
28
    using IndexSet = typename GridView::IndexSet;
    using Element  = typename GridView::template Codim<0>::Entity;
29
    using Estimates = std::vector<double>;
30
31
32
33
34
35

  public:
    /// Constructor.
    Marker() {}

    /// Constructor.
36
37
38
39
40
41
42
    Marker(std::string name, std::string component, Estimates const& est,
           std::shared_ptr<Grid> const& grid)
      : name_(name)
      , component_(component)
      , grid_(grid)
      , est_(est)
      , maximumMarking_(false)
43
    {
44
45
46
47
      Parameters::get(name_ + "->p", p_);
      Parameters::get(name_ + "->info", info_);
      Parameters::get(name_ + "->max refinement level", maxRefineLevel_);
      Parameters::get(name_ + "->min refinement level", minRefineLevel_);
48
49
50
51
52
    }

    /// destructor
    virtual ~Marker() {}

53
    /// Marks element with newMark. If \ref maximumMarking_ is set, the element
54
55
    /// is marked only if the new mark is bigger than the old one. The return
    /// value specifies whether the element has been marked, or not.
56
    void mark(Element const& elem, int newMark)
57
    {
58
      int oldMark = grid_->getMark(elem);
59

60
61
      if (!maximumMarking_ || (newMark > oldMark)) {
        bool marked = grid_->mark(newMark, elem);
62
63
        if (marked) {
          if (oldMark > 0) {
64
            elMarkRefine_--;
65
          } else if (oldMark < 0) {
66
            elMarkCoarsen_--;
67
68
69
          }

          if (newMark > 0) {
70
            elMarkRefine_++;
71
          } else if (newMark < 0) {
72
              elMarkCoarsen_++;
73
74
75
76
77
78
79
80
81
82
83
84
85
86
          }
        } else {
          msg("Marking failed");
        }
      }
    }

    /// Can be used by sub classes. Called before traversal.
    virtual void initMarking(AdaptInfo& adaptInfo);

    /// Can be used by sub classes. Called after traversal.
    virtual void finishMarking(AdaptInfo& adaptInfo);

    /// Marks one element.
87
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem);
88
89
90
91
92

    /// Marking of the mesh.
    virtual Flag markGrid(AdaptInfo& adaptInfo);

    /// Sets \ref maximumMarking.
93
    void setMaximumMarking(bool maxMark)
94
    {
95
      maximumMarking_ = maxMark;
96
97
    }

98
    int getElMarkRefine()
99
    {
100
      return elMarkRefine_;
101
102
    }

103
    int getElMarkCoarsen()
104
    {
105
      return elMarkCoarsen_;
106
107
108
    }

    /// Returns \ref name of the Marker
109
    std::string getName() const
110
    {
111
      return name_;
112
113
114
    }

    /// Creates a scalar marker depending on the strategy set in parameters.
115
    static std::shared_ptr<Marker<Traits> > createMarker(std::string name, std::string component, Estimates const& est, std::shared_ptr<Grid> const& grid);
116
117
118

  protected:
    /// Name of the scalar marker.
119
    std::string name_;
120

121
    /// Problem component to work on
122
    std::string component_;
123
124

    /// Pointer to the grid
125
    std::shared_ptr<Grid> grid_;
126
127

    /// Pointer to the local error estimates
128
    Estimates est_;
129
130

    /// estimation sum
131
    double estSum_ = 0.0;
132
133

    /// estmation maximum
134
    double estMax_ = 0.0;
135
136
137

    /// If true, elements are marked only if the new value is bigger than
    /// the current marking.
138
    bool maximumMarking_;
139
140
141

    /// Lower limit for error estimation, from which an element is marked for
    /// refinement
142
    double markRLimit_;
143
144
145

    /// Upper limit for error estimation, from which an element is marked for
    /// coarsening
146
    double markCLimit_;
147

148
149
    /// power in estimator norm TODO: is this info necessary in marker?
    double p_ = 2.0;
150
151

    /// Info level.
152
    int info_ = 10;
153
154

    /// Counter for elements marked for refinement
155
    int elMarkRefine_ = 0;
156
157

    /// Counter for elements marked for coarsening
158
    int elMarkCoarsen_ = 0;
159
160

    /// Maximal level of all elements.
161
    int maxRefineLevel_ = -1;
162
163

    /// Minimal level of all elements.
164
    int minRefineLevel_ = -1;
165

166
    bool refineAllowed_ = true;
167

168
    bool coarsenAllowed_ = false;
169
170
171
172
173
174
175
176
177
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Global refinement.
   */
178
  template <class Traits>
179
180
  class GRMarker
      : public Marker<Traits>
181
  {
182
183
184
185
    using Super = Marker<Traits>;
    using Grid      = typename Super::Grid;
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;
186
187
188

  public:
    /// Constructor.
189
190
191
    GRMarker(std::string name, std::string component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Marker<Traits>(std::move(name), std::move(component), est, grid)
192
193
194
    {}

    /// Implementation of Marker::markElement().
195
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem)
196
    {
197
      if (this->refineAllowed_)
198
199
200
201
202
203
204
205
206
207
208
        this->mark(elem, 1);
    }
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Maximum strategy.
   */
209

210
  template <class Traits>
211
212
  class MSMarker
      : public Marker<Traits>
213
  {
214
215
216
    using Super = Marker<Traits>;
    using Grid      = typename Super::Grid;
    using Estimates = typename Super::Estimates;
217
218
219

  public:
    /// Constructor.
220
221
222
    MSMarker(std::string name, std::string component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Super{std::move(name), std::move(component), est, grid}
223
    {
224
225
      Parameters::get(this->name_ + "->MSGamma", msGamma_);
      Parameters::get(this->name_ + "->MSGammaC", msGammaC_);
226
227
228
229
230
231
232
    }

    /// Implementation of Marker::initMarking().
    void initMarking(AdaptInfo& adaptInfo);

  protected:
    /// Marking parameter.
233
    double msGamma_ = 0.5;
234
235

    /// Marking parameter.
236
    double msGammaC_ = 0.1;
237
238
239
240
241
242
243
244
245
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Equidistribution strategy
   */
246

247
  template <class Traits>
248
249
  class ESMarker
      : public Marker<Traits>
250
  {
251
252
253
    using Super = Marker<Traits>;
    using Grid      = typename Super::Grid;
    using Estimates = typename Super::Estimates;
254
255
256

  public:
    /// Constructor.
257
258
259
    ESMarker(std::string name, std::string component, Estimates const& est,
             std::shared_ptr<Grid> const& grid)
      : Super{std::move(name), std::move(component), est, grid}
260
    {
261
262
      Parameters::get(this->name_ + "->ESTheta", esTheta_);
      Parameters::get(this->name_ + "->ESThetaC", esThetaC_);
263
264
265
266
267
268
269
    }

    /// Implementation of Marker::initMarking().
    virtual void initMarking(AdaptInfo& adaptInfo);

  protected:
    /// Marking parameter.
270
    double esTheta_ = 0.9;
271
272

    /// Marking parameter.
273
    double esThetaC_ = 0.2;
274
275
276
277
278
279
280
281
282
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Guaranteed error reduction strategy
   */
283

284
  template <class Traits>
285
286
  class GERSMarker
      : public Marker<Traits>
287
  {
288
289
290
291
    using Super = Marker<Traits>;
    using Grid      = typename Super::Grid;
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;
292
293
294

  public:
    /// Constructor.
295
296
297
    GERSMarker(std::string name, std::string component, Estimates const& est,
               std::shared_ptr<Grid> const& grid)
      : Super{std::move(name), std::move(component), est, grid}
298
    {
299
300
301
      Parameters::get(this->name_ + "->GERSThetaStar", gersThetaStar_);
      Parameters::get(this->name_ + "->GERSNu", gersNu_);
      Parameters::get(this->name_ + "->GERSThetaC", gersThetaC_);
302
303
304
305
306
307
308
    }

    /// Implementation of Marker::markGrid().
    virtual Flag markGrid(AdaptInfo& adaptInfo);

  protected:
    /// Refinement marking function.
309
    void markElementForRefinement(AdaptInfo& adaptInfo, Element const& elem);
310
311

    /// Coarsening marking function.
312
    void markElementForCoarsening(AdaptInfo& adaptInfo, Element const& elem);
313
314
315

  protected:
    /// Marking parameter.
316
    double gersSum_ = 0.0;
317
318

    /// Marking parameter.
319
    double oldErrSum_ = 0.0;
320
321

    /// Marking parameter.
322
    double gersThetaStar_ = 0.6;
323
324

    /// Marking parameter.
325
    double gersNu_ = 0.1;
326
327

    /// Marking parameter.
328
    double gersThetaC_ = 0.1;
329
330
331
  };

}
332

333
#include "Marker.inc.hpp"