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

Praetorius, Simon's avatar
Praetorius, Simon committed
3
#include <limits>
4
#include <string>
Praetorius, Simon's avatar
Praetorius, Simon committed
5
#include <utility>
6

Praetorius, Simon's avatar
Praetorius, Simon committed
7
#include <dune/common/std/optional.hh>
8
9
#include <dune/grid/common/grid.hh>

Praetorius, Simon's avatar
Praetorius, Simon committed
10
11
12
13
#include <amdis/common/ConceptsBase.hpp>

#include <amdis/gridfunctions/GridFunctionConcepts.hpp>

14
15
16
17
#include <amdis/AdaptInfo.hpp>
#include <amdis/Flag.hpp>
#include <amdis/Initfile.hpp>

Praetorius, Simon's avatar
Praetorius, Simon committed
18
19
namespace AMDiS
{
20
21
22
23
  /**
   * \ingroup Adaption
   *
   * \brief
Praetorius, Simon's avatar
Praetorius, Simon committed
24
   * Base class for all markers.
25
   */
Praetorius, Simon's avatar
Praetorius, Simon committed
26
  template <class Grid>
27
28
  class Marker
  {
29
  protected:
Praetorius, Simon's avatar
Praetorius, Simon committed
30
    using GridView = typename Grid::LeafGridView;
31
32
33
34
    using Element  = typename GridView::template Codim<0>::Entity;

  public:
    /// Constructor.
Praetorius, Simon's avatar
Praetorius, Simon committed
35
    Marker() = default;
36
37

    /// Constructor.
Praetorius, Simon's avatar
Praetorius, Simon committed
38
    Marker(std::string const& name, std::shared_ptr<Grid> const& grid)
39
40
      : name_(name)
      , grid_(grid)
41
    {
42
43
44
      Parameters::get(name + "->info", info_);
      Parameters::get(name + "->max refinement level", maxRefineLevel_);
      Parameters::get(name + "->min refinement level", minRefineLevel_);
45
46
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
47
48
    /// Destructor.
    virtual ~Marker() = default;
49

50
    /// Marks element with newMark. If \ref maximumMarking_ is set, the element
51
52
    /// 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.
Praetorius, Simon's avatar
Praetorius, Simon committed
53
    void mark(Element const& elem, int newMark);
54

Praetorius, Simon's avatar
Praetorius, Simon committed
55
    /// Called before traversal.
56
57
    virtual void initMarking(AdaptInfo& adaptInfo);

Praetorius, Simon's avatar
Praetorius, Simon committed
58
    /// Called after traversal.
59
60
61
    virtual void finishMarking(AdaptInfo& adaptInfo);

    /// Marks one element.
Praetorius, Simon's avatar
Praetorius, Simon committed
62
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem) = 0;
63

Praetorius, Simon's avatar
Praetorius, Simon committed
64
    /// Marking of the whole grid.
65
66
    virtual Flag markGrid(AdaptInfo& adaptInfo);

Praetorius, Simon's avatar
Praetorius, Simon committed
67
68
    /// Returns \ref elMarkRefine_ of the Marker
    int getElMarkRefine() const
69
    {
70
      return elMarkRefine_;
71
72
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
73
74
    /// Returns \ref elMarkCoarsen_ of the Marker
    int getElMarkCoarsen() const
75
    {
76
      return elMarkCoarsen_;
77
78
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
79
80
    /// Returns \ref name_ of the Marker
    std::string getName() const
81
    {
82
      return name_;
83
84
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
    /// Sets \ref maximumMarking_.
    void setMaximumMarking(bool maxMark)
    {
      maximumMarking_ = maxMark;
    }

    /// Sets \ref refineAllowed_.
    void allowRefinement(bool allow)
    {
      refineAllowed_ = allow;
    }

    /// Sets \ref coarsenAllowed_.
    void allowCoarsening(bool allow)
    {
      coarsenAllowed_ = allow;
    }
102
103

  protected:
Praetorius, Simon's avatar
Praetorius, Simon committed
104
    /// Name of the marker.
105
    std::string name_;
106
107

    /// Pointer to the grid
108
    std::shared_ptr<Grid> grid_;
109
110
111

    /// If true, elements are marked only if the new value is bigger than
    /// the current marking.
Praetorius, Simon's avatar
Praetorius, Simon committed
112
    bool maximumMarking_ = false;
113
114

    /// Info level.
Praetorius, Simon's avatar
Praetorius, Simon committed
115
    int info_ = 0;
116
117

    /// Counter for elements marked for refinement
118
    int elMarkRefine_ = 0;
119
120

    /// Counter for elements marked for coarsening
121
    int elMarkCoarsen_ = 0;
122
123

    /// Maximal level of all elements.
124
    int maxRefineLevel_ = std::numeric_limits<int>::max();
125
126

    /// Minimal level of all elements.
127
    int minRefineLevel_ = 0;
128

Praetorius, Simon's avatar
Praetorius, Simon committed
129
    /// Allow elements to be marked for refinement
130
    bool refineAllowed_ = true;
131

Praetorius, Simon's avatar
Praetorius, Simon committed
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
    /// Allow elements to be marked for coarsening
    bool coarsenAllowed_ = true;
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Base class for all estimator-based markers.
   */
  template <class Grid>
  class EstimatorMarker
      : public Marker<Grid>
  {
  protected:
    using Super = Marker<Grid>;
    using Element   = typename Super::Element;
    using Estimates = std::vector<double>;

  public:
    /// Constructor.
    EstimatorMarker() = default;

    /// Constructor.
    EstimatorMarker(std::string name, std::string component, Estimates const& est,
                    std::shared_ptr<Grid> const& grid)
      : Super{std::move(name), grid}
      , component_(std::move(component))
      , est_(est)
    {
      Parameters::get(this->name_ + "->p", p_);
    }

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

    /// Marks one element.
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem) override;

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

  protected:
    /// Problem component to work on
    std::string component_;

    /// Pointer to the local error estimates
    Estimates est_;

    /// Estimation sum
    double estSum_ = 0.0;

    /// Estmation maximum
    double estMax_ = 0.0;

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

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

    /// Power in estimator norm
    double p_ = 2.0;
199
200
201
202
203
204
205
206
207
  };


  /**
   * \ingroup Adaption
   *
   * \brief
   * Global refinement.
   */
Praetorius, Simon's avatar
Praetorius, Simon committed
208
  template <class Grid>
209
  class GRMarker
Praetorius, Simon's avatar
Praetorius, Simon committed
210
      : public EstimatorMarker<Grid>
211
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
212
    using Super = EstimatorMarker<Grid>;
213
214
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;
215
216
217

  public:
    /// Constructor.
Praetorius, Simon's avatar
Praetorius, Simon committed
218
    GRMarker(std::string const& name, std::string const& component, Estimates const& est,
219
             std::shared_ptr<Grid> const& grid)
Praetorius, Simon's avatar
Praetorius, Simon committed
220
      : Super{name, component, est, grid}
221
222
    {}

Praetorius, Simon's avatar
Praetorius, Simon committed
223
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem) override
224
    {
225
      if (this->refineAllowed_)
226
227
228
229
230
231
232
233
234
235
236
        this->mark(elem, 1);
    }
  };


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

Praetorius, Simon's avatar
Praetorius, Simon committed
238
  template <class Grid>
239
  class MSMarker
Praetorius, Simon's avatar
Praetorius, Simon committed
240
      : public EstimatorMarker<Grid>
241
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
242
    using Super = EstimatorMarker<Grid>;
243
    using Estimates = typename Super::Estimates;
244
245
246

  public:
    /// Constructor.
247
    MSMarker(std::string const& name, std::string component, Estimates const& est,
248
             std::shared_ptr<Grid> const& grid)
249
      : Super{name, std::move(component), est, grid}
250
    {
251
252
      Parameters::get(name + "->MSGamma", msGamma_);
      Parameters::get(name + "->MSGammaC", msGammaC_);
253
254
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
255
    virtual void initMarking(AdaptInfo& adaptInfo) override;
256
257
258

  protected:
    /// Marking parameter.
259
    double msGamma_ = 0.5;
260
261

    /// Marking parameter.
262
    double msGammaC_ = 0.1;
263
264
265
266
267
268
269
270
271
  };


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

Praetorius, Simon's avatar
Praetorius, Simon committed
273
  template <class Grid>
274
  class ESMarker
Praetorius, Simon's avatar
Praetorius, Simon committed
275
      : public EstimatorMarker<Grid>
276
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
277
    using Super = EstimatorMarker<Grid>;
278
    using Estimates = typename Super::Estimates;
279
280
281

  public:
    /// Constructor.
282
    ESMarker(std::string const& name, std::string component, Estimates const& est,
283
             std::shared_ptr<Grid> const& grid)
284
      : Super{name, std::move(component), est, grid}
285
    {
286
287
      Parameters::get(name + "->ESTheta", esTheta_);
      Parameters::get(name + "->ESThetaC", esThetaC_);
288
289
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
290
    virtual void initMarking(AdaptInfo& adaptInfo) override;
291
292
293

  protected:
    /// Marking parameter.
294
    double esTheta_ = 0.9;
295
296

    /// Marking parameter.
297
    double esThetaC_ = 0.2;
298
299
300
301
302
303
304
305
306
  };


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

Praetorius, Simon's avatar
Praetorius, Simon committed
308
  template <class Grid>
309
  class GERSMarker
Praetorius, Simon's avatar
Praetorius, Simon committed
310
      : public EstimatorMarker<Grid>
311
  {
Praetorius, Simon's avatar
Praetorius, Simon committed
312
    using Super = EstimatorMarker<Grid>;
313
314
    using Element   = typename Super::Element;
    using Estimates = typename Super::Estimates;
315
316
317

  public:
    /// Constructor.
318
    GERSMarker(std::string const& name, std::string component, Estimates const& est,
319
               std::shared_ptr<Grid> const& grid)
320
      : Super{name, std::move(component), est, grid}
321
    {
322
323
324
      Parameters::get(name + "->GERSThetaStar", gersThetaStar_);
      Parameters::get(name + "->GERSNu", gersNu_);
      Parameters::get(name + "->GERSThetaC", gersThetaC_);
325
326
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
327
    virtual Flag markGrid(AdaptInfo& adaptInfo) override;
328
329
330

  protected:
    /// Refinement marking function.
331
    void markElementForRefinement(AdaptInfo& adaptInfo, Element const& elem);
332
333

    /// Coarsening marking function.
334
    void markElementForCoarsening(AdaptInfo& adaptInfo, Element const& elem);
335
336
337

  protected:
    /// Marking parameter.
338
    double gersSum_ = 0.0;
339
340

    /// Marking parameter.
341
    double oldErrSum_ = 0.0;
342
343

    /// Marking parameter.
344
    double gersThetaStar_ = 0.6;
345
346

    /// Marking parameter.
347
    double gersNu_ = 0.1;
348
349

    /// Marking parameter.
350
    double gersThetaC_ = 0.1;
351
352
  };

Praetorius, Simon's avatar
Praetorius, Simon committed
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382

  /**
   * \ingroup Adaption
   *
   * \brief
   * Marker based on an indicator given as grid-function.
   *
   * On each element the grid-function is evaluated in the barycenter. The returned
   * values must be an integer (or convertible to an integer) indicating the
   * desired refinement level of this element. The element is marked for refinement
   * if the current level is < the desired level or coarsened, if >.
   *
   * \tparam Grid     An Implementation of the \ref Dune::Grid interface
   * \tparam GridFct  A grid-function with `Range` convertible to `int`
   */
  template <class Grid, class GridFct>
  class GridFunctionMarker
      : public Marker<Grid>
  {
    using Super = Marker<Grid>;
    using Element = typename Super::Element;

    template <class GF>
    using IsGridFunction = decltype(localFunction(std::declval<GF>()));

    static_assert(Dune::Std::is_detected<IsGridFunction,GridFct>::value, "");

  public:
    /// Constructor.
    template <class GF>
383
    GridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid, GF&& gf)
Praetorius, Simon's avatar
Praetorius, Simon committed
384
385
      : Super{name, grid}
      , gridFct_{makeGridFunction(std::forward<GF>(gf), grid->leafGridView())}
386
    {}
Praetorius, Simon's avatar
Praetorius, Simon committed
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405

    /// \brief Implementation of \ref Marker::markElement. Does nothing since marking is
    /// done in \ref markGrid().
    virtual void markElement(AdaptInfo& adaptInfo, Element const& elem) final {}

    /// Mark element for refinement, if grid-function \ref gridFct_ evaluates to
    /// a value larger than the current level and is marked for coarsening of
    /// the result is less than the current value.
    virtual Flag markGrid(AdaptInfo& adaptInfo) override;

  private:
    /// Indicator function for adaptation
    GridFct gridFct_;
  };

#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
  // Deduction guide for GridFunctionMarker class
  template <class Grid, class PreGridFct>
  GridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid,
406
                     PreGridFct&& preGridFct)
Praetorius, Simon's avatar
Praetorius, Simon committed
407
408
409
410
411
412
413
    -> GridFunctionMarker<Grid,
          std::decay_t<decltype(makeGridFunction(std::forward<PreGridFct>(preGridFct), grid->leafGridView()))>>;
#endif

  // Generator function for GridFunctionMarker class
  template <class Grid, class PreGridFct>
  auto makeGridFunctionMarker(std::string const& name, std::shared_ptr<Grid> const& grid,
414
                              PreGridFct&& preGridFct)
Praetorius, Simon's avatar
Praetorius, Simon committed
415
416
417
  {
    auto gridFct = makeGridFunction(std::forward<PreGridFct>(preGridFct), grid->leafGridView());
    using GridFct = decltype(gridFct);
418
    return GridFunctionMarker<Grid,GridFct>{name, grid, gridFct};
Praetorius, Simon's avatar
Praetorius, Simon committed
419
420
421
  }

} // end namespace AMDiS
422

423
#include "Marker.inc.hpp"