DataTransfer.hpp 4.26 KB
Newer Older
1
2
#pragma once

3
#include <map>
4
5
#include <memory>

6
7
#include <dune/grid/common/mcmgmapper.hh>

8
#include <amdis/Output.hpp>
9
#include <amdis/typetree/TreeContainer.hpp>
10
11
12

namespace AMDiS
{
13
14
15
16
17
  /**
   * \addtogroup Adaption
   * @{
   **/

18
  enum class DataTransferOperation {
19
20
    NO_OPERATION = 0,
    INTERPOLATE  = 1
21
  };
22
23


24
  /// \brief Interface for Containers allowing data transfer between grid changes.
25
26
27
28
29
30
31
32
33
34
35
  template <class Container>
  class DataTransferInterface
  {
  public:
    /// Virtual destructor
    virtual ~DataTransferInterface() = default;

    /// Collect data that is needed before grid adaption
    virtual void preAdapt(Container const& container, bool mightCoarsen) = 0;

    /// Interpolate data to new grid after grid adaption
36
    virtual void postAdapt(Container& container, bool refined) = 0;
37
38
39
  };


40
41
42
43
  /**
   * \brief Implementation of \ref DataTransferInterface that does not interpolate, but
   * just resizes the containers to the dimension of the basis
   **/
44
45
46
47
48
  template <class Container>
  class NoDataTransfer
      : public DataTransferInterface<Container>
  {
  public:
49
    void preAdapt(Container const& container, bool) override {}
50

51
    void postAdapt(Container& container, bool) override
52
    {
53
      container.resize();
54
55
56
57
    }
  };


58
59
60
61
62
63
64
65
  template <class Node, class Container, class Basis>
  class NodeDataTransfer;


  /** Data Transfer implementation for a single grid using interpolation
   *  Handles computations related to the geometric information of the grid and passes that to the
   *    underlying NodeDataTransfer classes
   */
66
  template <class Container, class Basis>
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
  class DataTransfer
      : public DataTransferInterface<Container>
  {
    using LocalView = typename Basis::LocalView;
    using Tree = typename LocalView::Tree;
    using GridView = typename Basis::GridView;
    using Grid = typename GridView::Grid;

    using Element = typename GridView::template Codim<0>::Entity;
    using Geometry = typename Element::Geometry;
    using LocalCoordinate = typename Geometry::LocalCoordinate;
    using IdType = typename Grid::LocalIdSet::IdType;

    template <class Node>
    using NodeElementData = typename NodeDataTransfer<Node, Container, Basis>::NodeElementData;
    using ElementData = TYPEOF(makeTreeContainer<NodeElementData>(std::declval<const Tree&>()));

  public:
    DataTransfer(std::shared_ptr<Basis> basis);

    /** Saves data contained in coeff in the PersistentContainer
     *  To be called after grid.preAdapt() and before grid.adapt()
     */
    void preAdapt(Container const& coeff, bool mightCoarsen) override;

    /** Unpacks data from the PersistentContainer
     *  To be called after grid.adapt() and before grid.postAdapt()
     */
95
96
    // [[expects: basis is update in gridView]]
    // [[expects: comm is updated in basis]]
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
    void postAdapt(Container& coeff, bool refined) override;

  private:
    /// The global basis
    std::shared_ptr<Basis const> basis_;

    /// Container with data that persists during grid adaptation
    using PersistentContainer = std::map<IdType, ElementData>;
    PersistentContainer persistentContainer_;

    /// Map leaf entities to unique index
    using Mapper = Dune::LeafMultipleCodimMultipleGeomTypeMapper<Grid>;
    Mapper mapper_;

    /// Data transfer on a single basis node
    template <class Node>
    using NDT = NodeDataTransfer<Node, Container, Basis>;
    using NodeDataTransferContainer = TYPEOF(makeTreeContainer<NDT>(std::declval<const Tree&>()));
    NodeDataTransferContainer nodeDataTransfer_;
  };

118
119
120
121
122
123
124
125
126

  /// Factory to create DataTransfer objects based on the \ref DataTransferOperation
  template <class Container>
  class DataTransferFactory
  {
    using Interface = DataTransferInterface<Container>;

  public:
    template <class Basis>
127
    static std::unique_ptr<Interface> create(DataTransferOperation op, std::shared_ptr<Basis> basis)
128
129
130
    {
      switch (op)
      {
131
        case DataTransferOperation::NO_OPERATION:
132
          return std::make_unique<NoDataTransfer<Container>>();
133
        case DataTransferOperation::INTERPOLATE:
134
          return std::make_unique<DataTransfer<Container, Basis>>(std::move(basis));
135
136
137
138
139
140
141
        default:
          error_exit("Invalid data transfer\n");
          return nullptr; // avoid warnings
      }
    }
  };

142
143
  /// @}

144
} // end namespace AMDiS
145
146

#include "DataTransfer.inc.hpp"