DOFVector.hpp 6.5 KB
Newer Older
1
2
3
4
5
6
7
8
9
#pragma once

#include <memory>
#include <utility>

#include <dune/common/shared_ptr.hh>

#include <amdis/DataTransfer.hpp>
#include <amdis/LinearAlgebra.hpp>
10
#include <amdis/Observer.hpp>
11
12
13
14
15
#include <amdis/common/Concepts.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/gridfunctions/GridFunction.hpp>
#include <amdis/typetree/TreePath.hpp>

16
17
18
19
20
21
22
23
24
namespace Dune
{
  namespace Functions
  {
    template <class PB>
    class DefaultGlobalBasis;
  }
}

25
26
namespace AMDiS
{
27
28
29
30
31
32
33
34
35
36
37
  // Forward declarations
  template <class PB>
  class ParallelGlobalBasis;

  namespace event
  {
    struct preAdapt;
    struct postAdapt;
    template <class B> struct basisUpdate;
  }

38
39
40
41
42
43
44
45
  /// \brief The basic container that stores a base vector and a corresponding basis
  /**
   * \tparam GB  Basis of the vector
   * \tparam T   Type of the coefficients
   **/
  template <class GB, class T = double>
  class DOFVector
      : public VectorBase<GB, VectorBackend<BackendTraits<GB,T>>>
46
      , public Observer<GB, event::preAdapt, event::basisUpdate<GB>, event::postAdapt>
47
48
49
  {
    using Self = DOFVector;
    using Super = VectorBase<GB, VectorBackend<BackendTraits<GB,T>>>;
50
    using Obs = Observer<GB, event::preAdapt, event::basisUpdate<GB>, event::postAdapt>;
51
52
53
54
55
56
57
58
59
60
61
62

  public:
    using Backend = VectorBackend<BackendTraits<GB,T>>;

    using GlobalBasis = GB;

    /// The index/size - type
    using size_type  = typename GlobalBasis::size_type;

    /// The type of the elements of the DOFVector
    using value_type = typename Backend::value_type;

63
64
    /// Wrapper for the implementation of the transfer of data during grid adaption
    using DataTransfer = DataTransferWrapper<Self>;
65
66
67

  public:
    /// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
68
69
70
71
72
73
    DOFVector(std::shared_ptr<GB> basis,
              DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : Super(basis)
      , Obs(basis)
      , dataTransfer_(op, basis)
      , basis_(basis)
74
75
    {}

76
77
78
79
80
81
    /// Constructor creating a new basis from a Dune::Functions::DefaultGlobalBasis.
    // TODO(FM): Replace explicit type with concept check
    DOFVector(Dune::Functions::DefaultGlobalBasis<typename GB::PreBasis>&& basis,
              DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(std::make_shared<GB>(std::move(basis)), op)
    {}
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

    template <class GB_,
      REQUIRES(Concepts::Similar<GB_,GB>)>
    DOFVector(GB_&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(Dune::wrap_or_move(FWD(basis)), op)
    {}


    template <class TreePath = RootTreePath>
    auto child(TreePath const& path = {})
    {
      auto&& tp = makeTreePath(path);
      return makeDOFVectorView(*this, tp);
    }

    template <class TreePath = RootTreePath>
    auto child(TreePath const& path = {}) const
    {
      auto&& tp = makeTreePath(path);
      return makeDiscreteFunction(*this, tp);
    }


    /// Interpolation of GridFunction to DOFVector, assuming that there is no
    /// reference to this DOFVector in the expression.
    /// See \ref DOFVectorView::interpolate_noalias
    template <class Expr, class Tag = tag::average>
    void interpolate_noalias(Expr&& expr, Tag strategy)
    {
      child().interpolate_noalias(FWD(expr), strategy);
    }

    /// Interpolation of GridFunction to DOFVector.
    /// See \ref DOFVectorView::interpolate
    template <class Expr, class Tag = tag::average>
    void interpolate(Expr&& expr, Tag strategy)
    {
      child().interpolate(FWD(expr), strategy);
    }

    /// Interpolation of GridFunction to DOFVector.
    /// See \ref DOFVectorView::interpolate
    template <class Expr>
    DOFVector& operator<<(Expr&& expr)
    {
      child().interpolate(FWD(expr));
      return *this;
    }


    /// Write DOFVector to file
    void backup(std::string const& filename);

    /// Read backup data from file
    void restore(std::string const& filename);


    /// Return the associated DataTransfer object
140
    DataTransfer const& dataTransfer() const
141
142
143
144
145
    {
      return dataTransfer_;
    }

    /// Return the associated DataTransfer object
146
    DataTransfer& dataTransfer()
147
148
149
150
151
    {
      return dataTransfer_;
    }

    /// Assign the DataTransfer object
152
    void setDataTransfer(DataTransfer const& dataTransfer)
153
154
155
156
    {
      dataTransfer_ = dataTransfer;
    }

157
158
159
160
    void setDataTransfer(DataTransfer&& dataTransfer)
    {
      dataTransfer_ = std::move(dataTransfer);
    }
161

162
163
    /// Create a new DataTransfer object based on the operation type
    void setDataTransfer(DataTransferOperation op)
164
    {
165
      setDataTransfer(newDataTransfer(op, basis_));
166
167
    }

168
169
170
    /// Override of Observer update(event::preAdapt) method. Redirects to preAdapt method of the
    /// \ref DataTransfer object.
    void update(event::preAdapt const& e) override
171
    {
172
      dataTransfer_.preAdapt(*this, e.mightCoarsen);
173
174
    }

175
176
177
    /// Override of Observer update(event::adapt) method. Redirects to adapt method of the
    /// \ref DataTransfer object.
    void update(event::basisUpdate<GB> const&) override
178
    {
179
180
      this->resize();
      dataTransfer_.adapt(*this);
181
182
    }

183
184
185
    /// Override of Observer update(event::postAdapt) method. Redirects to postAdapt method of the
    /// \ref DataTransfer object.
    void update(event::postAdapt const&) override
186
    {
187
      dataTransfer_.postAdapt(*this);
188
189
190
191
192
    }

  private:
    /// Data interpolation when the grid changes, set by default
    /// to \ref DataTransferOperation::INTERPOLATE.
193
194
195
    DataTransfer dataTransfer_;

    std::shared_ptr<GB> basis_;
196
197
198
199
200
  };

#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
  template <class GB>
  DOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
201
    -> DOFVector<ParallelGlobalBasis<typename Underlying_t<GB>::PreBasis>>;
202
203
204
205
206
207
208
209
210
211
212
213
#endif

  /// \brief Create a DOFVector from a basis.
  /**
   * This generator function accepts the basis as reference, temporary, or
   * shared_ptr. Internally the reference is wrapped into a non-destroying
   * shared_ptr and the temporary is moved into a new shared_ptr.
   *
   * The DataTransferOperation controls what is done during grid changes with the
   * DOFVector. The default is interpolation of the data to the new grid. See
   * \ref DataTransferOperation for more options.
   **/
214
215
216
  template <class ValueType = double, class GB>
  DOFVector<ParallelGlobalBasis<typename Underlying_t<GB>::PreBasis>, ValueType>
  makeDOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
217
218
219
220
221
222
223
  {
    return {FWD(basis), op};
  }

} // end namespace AMDiS

#include <amdis/DOFVector.inc.hpp>