DOFVector.hpp 6.3 KB
Newer Older
1
2
3
4
5
#pragma once

#include <memory>
#include <utility>

6
#include <dune/common/concept.hh>
7
#include <dune/common/shared_ptr.hh>
8
9
#include <dune/functions/functionspacebases/concepts.hh>

10
11
#include <amdis/DataTransfer.hpp>
#include <amdis/LinearAlgebra.hpp>
12
#include <amdis/Observer.hpp>
13
14
#include <amdis/common/Concepts.hpp>
#include <amdis/common/TypeTraits.hpp>
15
#include <amdis/functions/GlobalBasis.hpp>
16
17
18
19
20
21
22
23
#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
Praetorius, Simon's avatar
Praetorius, Simon committed
24
   * \tparam TraitsType  Collection of parameter for the linear-algebra backend
25
   **/
Praetorius, Simon's avatar
Praetorius, Simon committed
26
  template <class GB, class T = double, class TraitsType = BackendTraits<GB>>
27
  class DOFVector
Praetorius, Simon's avatar
Praetorius, Simon committed
28
      : public VectorFacade<T, TraitsType::template VectorImpl>
Praetorius, Simon's avatar
Praetorius, Simon committed
29
30
31
      , private Observer<event::preAdapt>
      , private Observer<event::adapt>
      , private Observer<event::postAdapt>
32
33
34
35
  {
    using Self = DOFVector;

  public:
36
    /// A global basis associated to the coefficients
37
38
    using GlobalBasis = GB;

39
40
41
    /// The internal coefficient vector storage
    using Coefficients = VectorFacade<T, TraitsType::template VectorImpl>;

42
43
44
45
    /// The index/size - type
    using size_type  = typename GlobalBasis::size_type;

    /// The type of the elements of the DOFVector
Praetorius, Simon's avatar
Praetorius, Simon committed
46
    using value_type = T;
47

48
49
    /// Wrapper for the implementation of the transfer of data during grid adaption
    using DataTransfer = DataTransferWrapper<Self>;
50

Praetorius, Simon's avatar
Praetorius, Simon committed
51
52
53
    /// Collection of parameter for the linear-algebra backend
    using Traits = TraitsType;

54
  public:
55
56
    /// (1) Constructor. Stores the shared_ptr of the basis and creates a new
    /// DataTransfer.
Praetorius, Simon's avatar
Praetorius, Simon committed
57
    DOFVector(std::shared_ptr<GB> const& basis,
58
              DataTransferOperation op = DataTransferOperation::INTERPOLATE)
Praetorius, Simon's avatar
Praetorius, Simon committed
59
      : Coefficients(*basis)
Praetorius, Simon's avatar
Praetorius, Simon committed
60
      , Observer<event::preAdapt>(basis->gridView().grid())
61
      , Observer<event::adapt>(*basis)
Praetorius, Simon's avatar
Praetorius, Simon committed
62
      , Observer<event::postAdapt>(basis->gridView().grid())
63
64
      , dataTransfer_(op, basis)
      , basis_(basis)
65
66
    {}

67
68
    /// (2) Constructor. Forwards to (1) by wrapping reference into non-destroying
    /// shared_ptr, see \ref Dune::wrap_or_move.
69
70
71
72
73
74
    template <class GB_,
      REQUIRES(Concepts::Similar<GB_,GB>)>
    DOFVector(GB_&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
      : DOFVector(Dune::wrap_or_move(FWD(basis)), op)
    {}

75
76
    /// (3) Constructor. Forwards to (1) by creating a new basis from a dune-functions
    /// pre-basis factory.
77
78
    template <class GV, class PBF,
      REQUIRES(Concepts::PreBasisFactory<PBF, GV, MultiIndex_t<PBF>>)>
79
80
    DOFVector(GV const& gridView, PBF const& preBasisFactory,
              DataTransferOperation op = DataTransferOperation::INTERPOLATE)
81
      : DOFVector(std::make_shared<GB>(gridView, preBasisFactory), op)
82
83
    {}

84
    /// Return the global basis
85
    std::shared_ptr<GlobalBasis const> const& basis() const
86
    {
87
      return basis_;
88
89
    }

90
    Coefficients const& coefficients() const
Praetorius, Simon's avatar
Praetorius, Simon committed
91
92
93
94
    {
      return static_cast<Coefficients const&>(*this);
    }

95
    Coefficients& coefficients()
Praetorius, Simon's avatar
Praetorius, Simon committed
96
97
98
99
    {
      return static_cast<Coefficients&>(*this);
    }

100
101
102
103
104
105
    /// Write DOFVector to file
    void backup(std::string const& filename);

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

106
107
108
109
110
111
112
113
114
115
116
117
118
    void resize()
    {
      Coefficients::resize(sizeInfo(*basis_));
    }

    using Coefficients::resize;

    void resizeZero()
    {
      Coefficients::resizeZero(sizeInfo(*basis_));
    }

    using Coefficients::resizeZero;
119
120

    /// Return the associated DataTransfer object
121
    DataTransfer const& dataTransfer() const
122
123
124
125
126
    {
      return dataTransfer_;
    }

    /// Return the associated DataTransfer object
127
    DataTransfer& dataTransfer()
128
129
130
131
132
    {
      return dataTransfer_;
    }

    /// Assign the DataTransfer object
133
    void setDataTransfer(DataTransfer const& dataTransfer)
134
135
136
137
    {
      dataTransfer_ = dataTransfer;
    }

138
139
140
141
    void setDataTransfer(DataTransfer&& dataTransfer)
    {
      dataTransfer_ = std::move(dataTransfer);
    }
142

143
144
    /// Create a new DataTransfer object based on the operation type
    void setDataTransfer(DataTransferOperation op)
145
    {
146
      setDataTransfer(DataTransfer(op, basis_));
147
148
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
149
150
  protected:

151
152
    /// Override of Observer update(event::preAdapt) method. Redirects to preAdapt
    /// method of the \ref DataTransfer object.
Praetorius, Simon's avatar
Praetorius, Simon committed
153
    void updateImpl(event::preAdapt e) override
154
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
155
      dataTransfer_.preAdapt(*this, e.value);
156
157
    }

158
159
    /// Override of Observer update(event::adapt) method. Redirects to adapt method
    /// of the \ref DataTransfer object.
Praetorius, Simon's avatar
Praetorius, Simon committed
160
    void updateImpl(event::adapt e) override
161
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
162
      assert(e.value);
163
      resize();
164
      dataTransfer_.adapt(*this);
165
166
    }

167
168
    /// Override of Observer update(event::postAdapt) method. Redirects to postAdapt
    /// method of the \ref DataTransfer object.
Praetorius, Simon's avatar
Praetorius, Simon committed
169
    void updateImpl(event::postAdapt) override
170
    {
171
      dataTransfer_.postAdapt(*this);
172
173
174
175
176
    }

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

Praetorius, Simon's avatar
Praetorius, Simon committed
179
    std::shared_ptr<GlobalBasis const> basis_;
180
181
  };

182
183

  // deduction guides
184
185
  template <class GB>
  DOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
186
    -> DOFVector<Underlying_t<GB>>;
187

188
  template <class GV, class PBF>
189
190
191
  DOFVector(GV const& gridView, PBF const& pbf,
            DataTransferOperation op = DataTransferOperation::INTERPOLATE)
    -> DOFVector<decltype(GlobalBasis{gridView,pbf})>;
192

193
194
195
196
197
198
199
200
201
202
203

  /// \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.
   **/
204
  template <class ValueType = double, class GB>
205
  DOFVector<Underlying_t<GB>, ValueType>
206
  makeDOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
207
208
209
210
211
212
213
  {
    return {FWD(basis), op};
  }

} // end namespace AMDiS

#include <amdis/DOFVector.inc.hpp>