AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
DOFVector.hpp
1 #pragma once
2 
3 #include <memory>
4 #include <utility>
5 
6 #include <dune/common/concept.hh>
7 #include <dune/common/shared_ptr.hh>
8 #include <dune/functions/functionspacebases/concepts.hh>
9 
10 #include <amdis/DataTransfer.hpp>
11 #include <amdis/LinearAlgebra.hpp>
12 #include <amdis/Observer.hpp>
13 #include <amdis/common/Concepts.hpp>
14 #include <amdis/common/TypeTraits.hpp>
15 #include <amdis/functions/GlobalBasis.hpp>
16 #include <amdis/typetree/TreePath.hpp>
17 
18 namespace AMDiS
19 {
21 
26  template <class GB, class T = double, class TraitsType = BackendTraits<GB>>
27  class DOFVector
28  : public VectorFacade<T, TraitsType::template VectorImpl>
29  , private Observer<event::preAdapt>
30  , private Observer<event::adapt>
31  , private Observer<event::postAdapt>
32  {
33  using Self = DOFVector;
34 
35  public:
37  using GlobalBasis = GB;
38 
41 
43  using size_type = typename GlobalBasis::size_type;
44 
46  using value_type = T;
47 
50 
52  using Traits = TraitsType;
53 
54  public:
57  DOFVector(std::shared_ptr<GB> const& basis,
58  DataTransferOperation op = DataTransferOperation::INTERPOLATE)
59  : Coefficients(*basis)
60  , Observer<event::preAdapt>(basis->gridView().grid())
61  , Observer<event::adapt>(*basis)
62  , Observer<event::postAdapt>(basis->gridView().grid())
63  , dataTransfer_(op, basis)
64  , basis_(basis)
65  {}
66 
69  template <class GB_,
70  REQUIRES(Concepts::Similar<GB_,GB>)>
71  DOFVector(GB_&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
72  : DOFVector(Dune::wrap_or_move(FWD(basis)), op)
73  {}
74 
77  template <class GV, class PBF,
78  REQUIRES(Concepts::PreBasisFactory<PBF, GV, MultiIndex_t<PBF>>)>
79  DOFVector(GV const& gridView, PBF const& preBasisFactory,
80  DataTransferOperation op = DataTransferOperation::INTERPOLATE)
81  : DOFVector(std::make_shared<GB>(gridView, preBasisFactory), op)
82  {}
83 
85  GlobalBasis const& basis() const
86  {
87  return *basis_;
88  }
89 
90  Coefficients const& coefficients() const
91  {
92  return static_cast<Coefficients const&>(*this);
93  }
94 
95  Coefficients& coefficients()
96  {
97  return static_cast<Coefficients&>(*this);
98  }
99 
101  void backup(std::string const& filename);
102 
104  void restore(std::string const& filename);
105 
106  void resize()
107  {
108  Coefficients::resize(sizeInfo(*basis_));
109  }
110 
111  using Coefficients::resize;
112 
113  void resizeZero()
114  {
115  Coefficients::resizeZero(sizeInfo(*basis_));
116  }
117 
119 
121  DataTransfer const& dataTransfer() const
122  {
123  return dataTransfer_;
124  }
125 
128  {
129  return dataTransfer_;
130  }
131 
134  {
135  dataTransfer_ = dataTransfer;
136  }
137 
139  {
140  dataTransfer_ = std::move(dataTransfer);
141  }
142 
144  void setDataTransfer(DataTransferOperation op)
145  {
146  setDataTransfer(DataTransfer(op, basis_));
147  }
148 
149  protected:
150 
153  void updateImpl(event::preAdapt e) override
154  {
155  dataTransfer_.preAdapt(*this, e.value);
156  }
157 
160  void updateImpl(event::adapt e) override
161  {
162  assert(e.value);
163  resize();
164  dataTransfer_.adapt(*this);
165  }
166 
170  {
171  dataTransfer_.postAdapt(*this);
172  }
173 
174  private:
177  DataTransfer dataTransfer_;
178 
179  std::shared_ptr<GlobalBasis const> basis_;
180  };
181 
182 
183  // deduction guides
184  template <class GB>
185  DOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
187 
188  template <class GV, class PBF>
189  DOFVector(GV const& gridView, PBF const& pbf,
190  DataTransferOperation op = DataTransferOperation::INTERPOLATE)
191  -> DOFVector<decltype(GlobalBasis{gridView,pbf})>;
192 
193 
195 
204  template <class ValueType = double, class GB>
205  DOFVector<Underlying_t<GB>, ValueType>
206  makeDOFVector(GB&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
207  {
208  return {FWD(basis), op};
209  }
210 
211 } // end namespace AMDiS
212 
213 #include <amdis/DOFVector.inc.hpp>
GB GlobalBasis
A global basis associated to the coefficients.
Definition: DOFVector.hpp:37
DOFVector(GB_ &&basis, DataTransferOperation op=DataTransferOperation::INTERPOLATE)
Definition: DOFVector.hpp:71
T value_type
The type of the elements of the DOFVector.
Definition: DOFVector.hpp:46
The basic container that stores a base vector and a corresponding basis.
Definition: LinearSolverInterface.hpp:23
DataTransferWrapper< Self > DataTransfer
Wrapper for the implementation of the transfer of data during grid adaption.
Definition: DOFVector.hpp:49
Definition: AdaptiveGrid.hpp:373
Definition: FieldMatVec.hpp:12
DOFVector(std::shared_ptr< GB > const &basis, DataTransferOperation op=DataTransferOperation::INTERPOLATE)
Definition: DOFVector.hpp:57
void updateImpl(event::postAdapt) override
Definition: DOFVector.hpp:169
void updateImpl(event::adapt e) override
Definition: DOFVector.hpp:160
The basic container that stores a base vector and a corresponding basis.
Definition: DOFVector.hpp:27
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
DataTransfer & dataTransfer()
Return the associated DataTransfer object.
Definition: DOFVector.hpp:127
Implementation of the ObserverInterface.
Definition: Observer.hpp:100
DOFVector< Underlying_t< GB >, ValueType > makeDOFVector(GB &&basis, DataTransferOperation op=DataTransferOperation::INTERPOLATE)
Create a DOFVector from a basis.
Definition: DOFVector.hpp:206
void setDataTransfer(DataTransfer const &dataTransfer)
Assign the DataTransfer object.
Definition: DOFVector.hpp:133
void resize(SizeInfo const &sizeInfo)
Resize the vector to the size of the basis.
Definition: VectorFacade.hpp:101
Definition: Observer.hpp:25
Definition: Observer.hpp:30
void resizeZero(SizeInfo const &sizeInfo)
Resize the vector to the size of the basis and set to zero.
Definition: VectorFacade.hpp:108
void setDataTransfer(DataTransferOperation op)
Create a new DataTransfer object based on the operation type.
Definition: DOFVector.hpp:144
GlobalBasis const & basis() const
Return the global basis.
Definition: DOFVector.hpp:85
Definition: Observer.hpp:19
DataTransfer const & dataTransfer() const
Return the associated DataTransfer object.
Definition: DOFVector.hpp:121
void backup(std::string const &filename)
Write DOFVector to file.
Definition: DOFVector.inc.hpp:15
DOFVector(GV const &gridView, PBF const &preBasisFactory, DataTransferOperation op=DataTransferOperation::INTERPOLATE)
Definition: DOFVector.hpp:79
typename GlobalBasis::size_type size_type
The index/size - type.
Definition: DOFVector.hpp:43
void updateImpl(event::preAdapt e) override
Definition: DOFVector.hpp:153
void restore(std::string const &filename)
Read backup data from file.
Definition: DOFVector.inc.hpp:40