dune-vtk 2.8
Loading...
Searching...
No Matches
lagrangegridfunction.hh
Go to the documentation of this file.
1#pragma once
2
3#include <optional>
4#include <vector>
5
6#include <dune/common/dynvector.hh>
7#include <dune/common/version.hh>
8#include <dune/localfunctions/lagrange.hh>
13
14namespace Dune
15{
16 namespace Vtk
17 {
18 template <class GridType>
20
23 template <class GridType, class FieldType, class Context>
25 {
26 using Grid = GridType;
27 using Field = FieldType;
28
29 using Factory = GridFactory<Grid>;
31
32 public:
33 struct EntitySet
34 {
35 using Grid = GridType;
36 using Element = typename GridType::template Codim<0>::Entity;
37 using LocalCoordinate = typename Element::Geometry::LocalCoordinate;
38 using GlobalCoordinate = typename Element::Geometry::GlobalCoordinate;
39 };
40
42 using Range = DynamicVector<Field>;
44
45 private:
46 template <class LC>
47 class PointDataLocalFunction
48 {
49 using LFE = LagrangeLocalFiniteElement<LagrangePointSet, LC::dimension, FieldType, FieldType>;
50 using LB = typename LFE::Traits::LocalBasisType;
51
52 public:
53 using LocalContext = LC;
54 using Domain = typename LC::Geometry::LocalCoordinate;
55 using Range = DynamicVector<Field>;
56 using Signature = Range(Domain);
57
58 public:
60 PointDataLocalFunction (GridCreator const* creator, std::vector<Field> const* values, unsigned int comp,
61 std::vector<std::uint8_t> const* types,
62 std::vector<std::int64_t> const* offsets,
63 std::vector<std::int64_t> const* connectivity)
64 : creator_(creator)
65 , values_(values)
66 , comp_(comp)
67 , types_(types)
68 , offsets_(offsets)
69 , connectivity_(connectivity)
70 {}
71
72 PointDataLocalFunction () = default;
73
75
81 void bind (LocalContext const& element)
82 {
83 unsigned int insertionIndex = creator_->factory().insertionIndex(element);
84
85 std::int64_t shift = (insertionIndex == 0 ? 0 : (*offsets_)[insertionIndex-1]);
86 std::int64_t numNodes = (*offsets_)[insertionIndex] - shift;
87 [[maybe_unused]] std::int64_t maxNumNodes = numLagrangePoints(element.type(), 20);
88 VTK_ASSERT(numNodes > 0 && numNodes < maxNumNodes);
89
90 int order = creator_->order(element.type(), numNodes);
91 VTK_ASSERT(order > 0 && order < 20);
92
93 // construct a local finite-element with the corresponding order and Lagrange points
94 // as stored in the file
95 lfe_.emplace(LFE{element.type(), (unsigned int)(order)});
96 lgeo_.emplace(creator_->localGeometry(element));
97
98 // collect values on lagrange nodes
99 localValues_.resize(numNodes);
100 for (std::int64_t i = shift, i0 = 0; i < (*offsets_)[insertionIndex]; ++i, ++i0) {
101 std::int64_t idx = (*connectivity_)[i];
102 DynamicVector<Field>& v = localValues_[i0];
103 v.resize(comp_);
104 for (unsigned int j = 0; j < comp_; ++j)
105 v[j] = (*values_)[comp_*idx + j];
106 }
107 }
108
110 void unbind ()
111 {
112 lfe_.reset();
113 lgeo_.reset();
114 }
115
118 // NOTE: do we need to transform the local coordinates?
119 Range operator() (Domain const& local) const
120 {
121 assert(!!lfe_);
122 auto const& lb = lfe_->localBasis();
123 lb.evaluateFunction(lgeo_->global(local), shapeValues_);
124 assert(shapeValues_.size() == localValues_.size());
125
126 Range y(comp_, Field(0));
127 for (std::size_t i = 0; i < shapeValues_.size(); ++i)
128 y.axpy(shapeValues_[i], localValues_[i]);
129
130 return y;
131 }
132
133 private:
134 GridCreator const* creator_ = nullptr;
135 std::vector<Field> const* values_ = nullptr;
136 unsigned int comp_;
137 std::vector<std::uint8_t> const* types_ = nullptr;
138 std::vector<std::int64_t> const* offsets_ = nullptr;
139 std::vector<std::int64_t> const* connectivity_ = nullptr;
140
141 // Local Finite-Element
142 std::optional<LFE> lfe_ = std::nullopt;
143 std::optional<typename GridCreator::LocalGeometry> lgeo_ = std::nullopt;
144
145 // cache of local values
146 std::vector<DynamicVector<Field>> localValues_;
147 mutable std::vector<typename LB::Traits::RangeType> shapeValues_;
148 };
149
151 template <class LC>
152 class CellDataLocalFunction
153 {
154 public:
155 using LocalContext = LC;
156 using Domain = typename LC::Geometry::LocalCoordinate;
157 using Range = DynamicVector<Field>;
158 using Signature = Range(Domain);
159
160 public:
162 CellDataLocalFunction (GridCreator const* creator, std::vector<Field> const* values, unsigned int comp,
163 std::vector<std::uint8_t> const* /*types*/,
164 std::vector<std::int64_t> const* /*offsets*/,
165 std::vector<std::int64_t> const* /*connectivity*/)
166 : creator_(creator)
167 , values_(values)
168 , comp_(comp)
169 {}
170
171 CellDataLocalFunction () = default;
172
175 void bind (LocalContext const& element)
176 {
177 unsigned int idx = creator_->factory().insertionIndex(element);
178
179 // collect values on cells
180 DynamicVector<Field>& v = localValue_;
181 v.resize(comp_);
182
183 for (unsigned int j = 0; j < comp_; ++j)
184 v[j] = (*values_)[comp_*idx + j];
185 }
186
188 void unbind ()
189 {}
190
193 Range operator() (Domain const& local) const
194 {
195 return localValue_;
196 }
197
198 private:
199 GridCreator const* creator_ = nullptr;
200 std::vector<Field> const* values_ = nullptr;
201 unsigned int comp_;
202
203 // cache of local values
204 DynamicVector<Field> localValue_;
205 };
206
208 template <class LC>
209 using LocalFunction = std::conditional_t< std::is_same_v<Context,PointContext>,
210 PointDataLocalFunction<LC>,
211 CellDataLocalFunction<LC>>;
212
213 public:
216 LagrangeGridFunction (GridCreator const& creator, std::vector<Field> const& values,
217 std::string name ,unsigned int ncomps, Vtk::DataTypes dataType,
218 std::vector<std::uint8_t> const& types,
219 std::vector<std::int64_t> const& offsets,
220 std::vector<std::int64_t> const& connectivity)
221 : creator_(&creator)
222 , values_(&values)
223 , name_(std::move(name))
224 , ncomps_(ncomps)
225 , dataType_(dataType)
226 , types_(&types)
227 , offsets_(&offsets)
228 , connectivity_(&connectivity)
229 {}
230
232
234 Range operator() (Domain const& global) const
235 {
236 DUNE_THROW(Dune::NotImplemented, "Evaluation in global coordinates not implemented.");
237 return Range(ncomps_, 0);
238 }
239
241 EntitySet const& entitySet () const
242 {
243 return entitySet_;
244 }
245
246 std::string const& name () const
247 {
248 return name_;
249 }
250
251 int numComponents () const
252 {
253 return ncomps_;
254 }
255
257 {
258 return dataType_;
259 }
260
263 friend LocalFunction<typename EntitySet::Element> localFunction (LagrangeGridFunction const& gf)
264 {
265 return {gf.creator_, gf.values_, gf.ncomps_, gf.types_, gf.offsets_, gf.connectivity_};
266 }
267
268 private:
269 GridCreator const* creator_ = nullptr;
270 std::vector<Field> const* values_ = nullptr;
271 std::string name_ = "GridFunction";
272 unsigned int ncomps_ = 0;
274 std::vector<std::uint8_t> const* types_ = nullptr;
275 std::vector<std::int64_t> const* offsets_ = nullptr;
276 std::vector<std::int64_t> const* connectivity_ = nullptr;
277
278 EntitySet entitySet_;
279 };
280
281 } // end namespace Vtk
282} // end namespace Dune
Macro for wrapping error checks and throwing exceptions.
#define VTK_ASSERT(cond)
check if condition cond holds; otherwise, throw a VtkError.
Definition: errors.hh:29
Definition: writer.hh:13
LagrangeGridCreator(GridFactory< Grid > &) -> LagrangeGridCreator< Grid >
DataTypes
Definition: types.hh:52
Definition: lagrangegridcreator.hh:40
Grid-function representing values from a VTK file with local Lagrange interpolation of the values sto...
Definition: lagrangegridfunction.hh:25
DynamicVector< Field > Range
Definition: lagrangegridfunction.hh:42
EntitySet const & entitySet() const
Return a type that defines the element that can be iterated.
Definition: lagrangegridfunction.hh:241
Range(Domain) Signature
Definition: lagrangegridfunction.hh:43
int numComponents() const
Definition: lagrangegridfunction.hh:251
friend LocalFunction< typename EntitySet::Element > localFunction(LagrangeGridFunction const &gf)
Definition: lagrangegridfunction.hh:263
std::string const & name() const
Definition: lagrangegridfunction.hh:246
typename EntitySet::GlobalCoordinate Domain
Definition: lagrangegridfunction.hh:41
Range operator()(Domain const &global) const
Global evaluation. Not supported!
Definition: lagrangegridfunction.hh:234
LagrangeGridFunction(GridCreator const &creator, std::vector< Field > const &values, std::string name, unsigned int ncomps, Vtk::DataTypes dataType, std::vector< std::uint8_t > const &types, std::vector< std::int64_t > const &offsets, std::vector< std::int64_t > const &connectivity)
Definition: lagrangegridfunction.hh:216
Vtk::DataTypes dataType() const
Definition: lagrangegridfunction.hh:256
Definition: lagrangegridfunction.hh:34
typename Element::Geometry::GlobalCoordinate GlobalCoordinate
Definition: lagrangegridfunction.hh:38
typename Element::Geometry::LocalCoordinate LocalCoordinate
Definition: lagrangegridfunction.hh:37
GridType Grid
Definition: lagrangegridfunction.hh:35
typename GridType::template Codim< 0 >::Entity Element
Definition: lagrangegridfunction.hh:36