DOFVector.hpp 8.33 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
28
29
30
31
32
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
78
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
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
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
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
#pragma once

#include <memory>
#include <utility>

#include <dune/common/shared_ptr.hh>

#include <amdis/DataTransfer.hpp>
#include <amdis/DOFVectorInterface.hpp>
#include <amdis/GridTransferManager.hpp>
#include <amdis/LinearAlgebra.hpp>
#include <amdis/common/Concepts.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/gridfunctions/GridFunction.hpp>
#include <amdis/typetree/TreePath.hpp>

namespace AMDiS
{
  /// \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>>>
      , public DOFVectorInterface
  {
    using Self = DOFVector;
    using Super = VectorBase<GB, VectorBackend<BackendTraits<GB,T>>>;

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

    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;

    /// Defines an interface to transfer the data during grid adaption
    using DataTransfer = DataTransferInterface<Self>;

    /// A creator for a concrete data transfer object, depending on \ref DataTransferOperation
    using DataTransferFactory = AMDiS::DataTransferFactory<Self>;

  public:
    /// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
    DOFVector(std::shared_ptr<GlobalBasis> basis, std::shared_ptr<Comm> comm, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : Super(basis, std::move(comm))
      , dataTransfer_(DataTransferFactory::create(op, basis))
    {
      attachToGridTransfer();
    }

    /// Forwarding constructor that wraps the arguments into shared_ptr
    template <class GB_, class Comm_,
      REQUIRES(Concepts::Similar<GB_,GB>),
      REQUIRES(Concepts::Similar<Comm_,Comm>)>
    DOFVector(GB_&& basis, Comm_&& comm, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(Dune::wrap_or_move(FWD(basis)), Dune::wrap_or_move(FWD(comm)), op)
    {}


    /// Construct the DOFVector from a basis by first creating a comm object. This constructor
    /// Registers the basis and comm object into the \ref GridTransferManager so that it gets updated
    /// on grid changes.
    DOFVector(std::shared_ptr<GlobalBasis> basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(basis, CommunicationCreator<Comm>::create(*basis), op)
    {
      GridTransferManager::attach(this->basis()->gridView().grid(), *this->comm(), [c=this->comm(), b=this->basis()]() -> void {
        b->update(b->gridView());
        c->update(*b);
      });
    }

    /// Forwarding constructor that wraps the arguments into shared_ptr
    template <class GB_,
      REQUIRES(Concepts::Similar<GB_,GB>)>
    DOFVector(GB_&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(Dune::wrap_or_move(FWD(basis)), op)
    {}

    /// Copy constructor
    DOFVector(Self const& that)
      : Super(static_cast<Super const&>(that))
      , dataTransfer_(that.dataTransfer_)
    {
      attachToGridTransfer();
    }

    /// Move constructor
    DOFVector(Self&& that)
      : Super(static_cast<Super&&>(that))
      , dataTransfer_(std::move(that.dataTransfer_))
    {
      attachToGridTransfer();
    }

    /// Destructor
    ~DOFVector() override
    {
      detachFromGridTransfer();
    }

    /// Copy assignment operator
    Self& operator=(Self const& that)
    {
      detachFromGridTransfer();
      Super::operator=(that);
      dataTransfer_ = that.dataTransfer_;
      attachToGridTransfer();
      return *this;
    }

    /// Move assignment
    Self& operator=(Self&& that) = default;

    void resize() override
    {
      Super::resize();
    }


    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
    std::shared_ptr<DataTransfer const> dataTransfer() const
    {
      return dataTransfer_;
    }

    /// Return the associated DataTransfer object
    std::shared_ptr<DataTransfer> dataTransfer()
    {
      return dataTransfer_;
    }

    /// Create a new DataTransfer object based on the operation type
    void setDataTransfer(DataTransferOperation op)
    {
      dataTransfer_ = DataTransferFactory::create(op, this->basis());
    }

    /// Assign the DataTransfer object
    void setDataTransfer(std::shared_ptr<DataTransfer> const& dataTransfer)
    {
      dataTransfer_ = dataTransfer;
    }


    /// Implementation of \ref DOFVectorInterface::preAdapt
    /// Redirects to a \ref DataTransfer object.
    void preAdapt(bool mightCoarsen) override
    {
      dataTransfer_->preAdapt(*this, mightCoarsen);
    }

    /// Implementation of \ref DOFVectorInterface::postAdapt
    /// Redirects to a \ref DataTransfer object.
    void postAdapt(bool refined) override
    {
      dataTransfer_->postAdapt(*this, refined);
    }

  private:
    // register this DOFVector and its basis to the DataTransfer
    void attachToGridTransfer()
    {
      GridTransferManager::attach(this->basis()->gridView().grid(), *this);
    }

    // deregister this DOFVector and its basis from the DataTransfer
    void detachFromGridTransfer()
    {
      GridTransferManager::detach(this->basis()->gridView().grid(), *this);
    }

  private:
    /// Data interpolation when the grid changes, set by default
    /// to \ref DataTransferOperation::INTERPOLATE.
    std::shared_ptr<DataTransfer> dataTransfer_;
  };

#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
  template <class GB, class C>
  DOFVector(GB&& basis, C&& comm, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
    -> DOFVector<Underlying_t<GB>>;

  template <class GB>
  DOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
    -> DOFVector<Underlying_t<GB>>;
#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.
   **/
  template <class ValueType = double, class GlobalBasis, class Communication>
  DOFVector<Underlying_t<GlobalBasis>, ValueType>
  makeDOFVector(GlobalBasis&& basis, Communication&& comm, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
  {
    return {FWD(basis), FWD(comm), op};
  }

  template <class ValueType = double, class GlobalBasis>
  DOFVector<Underlying_t<GlobalBasis>, ValueType>
  makeDOFVector(GlobalBasis&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
  {
    return {FWD(basis), op};
  }


} // end namespace AMDiS

#include <amdis/DOFVector.inc.hpp>