6 #include <dune/common/dynvector.hh>
7 #include <dune/localfunctions/lagrange.hh>
17 template <
class Gr
idType>
22 template <
class Gr
idType,
class FieldType,
class Context>
25 using Grid = GridType;
26 using Field = FieldType;
28 using Factory = GridFactory<Grid>;
35 using Element =
typename GridType::template Codim<0>::Entity;
41 using Range = DynamicVector<Field>;
46 class PointDataLocalFunction
48 using LFE = LagrangeLocalFiniteElement<LagrangePointSet, LC::dimension, FieldType, FieldType>;
49 using LB =
typename LFE::Traits::LocalBasisType;
52 using LocalContext = LC;
53 using Domain =
typename LC::Geometry::LocalCoordinate;
54 using Range = DynamicVector<Field>;
59 PointDataLocalFunction (
GridCreator const& creator, std::vector<Field>
const& values,
unsigned int comp,
60 std::vector<std::uint8_t>
const& types,
61 std::vector<std::int64_t>
const& offsets,
62 std::vector<std::int64_t>
const& connectivity)
68 , connectivity_(connectivity)
78 void bind (LocalContext
const& element)
80 unsigned int insertionIndex = creator_.factory().insertionIndex(element);
82 std::int64_t shift = (insertionIndex == 0 ? 0 : offsets_[insertionIndex-1]);
83 std::int64_t numNodes = offsets_[insertionIndex] - shift;
84 [[maybe_unused]] std::int64_t maxNumNodes = numLagrangePoints(element.type().id(), element.type().dim(), 20);
85 VTK_ASSERT(numNodes > 0 && numNodes < maxNumNodes);
87 int order = creator_.order(element.type(), numNodes);
92 lfe_.emplace(LFE{element.type(), (
unsigned int)(order)});
93 lgeo_.emplace(creator_.localGeometry(element));
96 localValues_.resize(numNodes);
97 for (std::int64_t i = shift, i0 = 0; i < offsets_[insertionIndex]; ++i, ++i0) {
98 std::int64_t idx = connectivity_[i];
99 DynamicVector<Field>& v = localValues_[i0];
101 for (
unsigned int j = 0; j < comp_; ++j)
102 v[j] = values_[comp_*idx + j];
119 auto const& lb = lfe_->localBasis();
120 lb.evaluateFunction(lgeo_->global(local), shapeValues_);
121 assert(shapeValues_.size() == localValues_.size());
123 Range y(comp_, Field(0));
124 for (std::size_t i = 0; i < shapeValues_.size(); ++i)
125 y.axpy(shapeValues_[i], localValues_[i]);
131 GridCreator
const& creator_;
132 std::vector<Field>
const& values_;
134 std::vector<std::uint8_t>
const& types_;
135 std::vector<std::int64_t>
const& offsets_;
136 std::vector<std::int64_t>
const& connectivity_;
139 std::optional<LFE> lfe_ = std::nullopt;
140 std::optional<typename GridCreator::LocalGeometry> lgeo_ = std::nullopt;
143 std::vector<DynamicVector<Field>> localValues_;
144 mutable std::vector<typename LB::Traits::RangeType> shapeValues_;
149 class CellDataLocalFunction
152 using LocalContext = LC;
153 using Domain =
typename LC::Geometry::LocalCoordinate;
154 using Range = DynamicVector<Field>;
159 CellDataLocalFunction (GridCreator
const& creator, std::vector<Field>
const& values,
unsigned int comp,
160 std::vector<std::uint8_t>
const& ,
161 std::vector<std::int64_t>
const& ,
162 std::vector<std::int64_t>
const& )
163 : factory_(creator.factory())
170 void bind (LocalContext
const& element)
172 unsigned int idx = factory_.insertionIndex(element);
175 DynamicVector<Field>& v = localValue_;
178 for (
unsigned int j = 0; j < comp_; ++j)
179 v[j] = values_[comp_*idx + j];
194 Factory
const& factory_;
195 std::vector<Field>
const& values_;
199 DynamicVector<Field> localValue_;
204 using LocalFunction = std::conditional_t< std::is_same_v<Context,PointContext>,
205 PointDataLocalFunction<LC>,
206 CellDataLocalFunction<LC>>;
212 std::vector<std::uint8_t>
const& types,
213 std::vector<std::int64_t>
const& offsets,
214 std::vector<std::int64_t>
const& connectivity)
220 , connectivity_(connectivity)
226 DUNE_THROW(Dune::NotImplemented,
"Evaluation in global coordinates not implemented.");
227 return Range(comp_, 0);
240 return {gf.creator_, gf.values_, gf.comp_, gf.types_, gf.offsets_, gf.connectivity_};
244 GridCreator
const& creator_;
245 std::vector<Field>
const& values_;
247 std::vector<std::uint8_t>
const& types_;
248 std::vector<std::int64_t>
const& offsets_;
249 std::vector<std::int64_t>
const& connectivity_;
251 EntitySet entitySet_;
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: datacollectorinterface.hh:9
LagrangeGridCreator(GridFactory< Grid > &) -> LagrangeGridCreator< Grid >
Definition: lagrangegridcreator.hh:40
Grid-function representing values from a VTK file with local Lagrange interpolation of the values sto...
Definition: lagrangegridfunction.hh:24
DynamicVector< Field > Range
Definition: lagrangegridfunction.hh:41
friend LocalFunction< typename EntitySet::Element > localFunction(LagrangeGridFunction const &gf)
Definition: lagrangegridfunction.hh:238
LagrangeGridFunction(GridCreator const &creator, std::vector< Field > const &values, unsigned int comp, std::vector< std::uint8_t > const &types, std::vector< std::int64_t > const &offsets, std::vector< std::int64_t > const &connectivity)
Definition: lagrangegridfunction.hh:211
Range(Domain) Signature
Definition: lagrangegridfunction.hh:42
typename EntitySet::GlobalCoordinate Domain
Definition: lagrangegridfunction.hh:40
Range operator()(Domain const &global) const
Global evaluation. Not supported!
Definition: lagrangegridfunction.hh:224
EntitySet const & entitySet() const
Return a type that defines the element that can be iterated.
Definition: lagrangegridfunction.hh:231
Definition: lagrangegridfunction.hh:33
typename Element::Geometry::GlobalCoordinate GlobalCoordinate
Definition: lagrangegridfunction.hh:37
typename Element::Geometry::LocalCoordinate LocalCoordinate
Definition: lagrangegridfunction.hh:36
GridType Grid
Definition: lagrangegridfunction.hh:34
typename GridType::template Codim< 0 >::Entity Element
Definition: lagrangegridfunction.hh:35