6 #include <dune/common/dynvector.hh>
7 #include <dune/common/version.hh>
8 #include <dune/localfunctions/lagrange.hh>
18 template <
class Gr
idType>
23 template <
class Gr
idType,
class FieldType,
class Context>
26 using Grid = GridType;
27 using Field = FieldType;
29 using Factory = GridFactory<Grid>;
36 using Element =
typename GridType::template Codim<0>::Entity;
42 using Range = DynamicVector<Field>;
47 class PointDataLocalFunction
49 using LFE = LagrangeLocalFiniteElement<LagrangePointSet, LC::dimension, FieldType, FieldType>;
50 using LB =
typename LFE::Traits::LocalBasisType;
53 using LocalContext = LC;
54 using Domain =
typename LC::Geometry::LocalCoordinate;
55 using Range = DynamicVector<Field>;
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)
69 , connectivity_(connectivity)
72 PointDataLocalFunction () =
default;
81 void bind (LocalContext
const& element)
83 unsigned int insertionIndex = creator_->factory().insertionIndex(element);
85 std::int64_t shift = (insertionIndex == 0 ? 0 : (*offsets_)[insertionIndex-1]);
86 std::int64_t numNodes = (*offsets_)[insertionIndex] - shift;
87 #if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
88 [[maybe_unused]] std::int64_t maxNumNodes = numLagrangePoints(element.type().id(), element.type().dim(), 20);
90 [[maybe_unused]] std::int64_t maxNumNodes = numLagrangePoints(element.type(), 20);
92 VTK_ASSERT(numNodes > 0 && numNodes < maxNumNodes);
94 int order = creator_->order(element.type(), numNodes);
99 lfe_.emplace(LFE{element.type(), (
unsigned int)(order)});
100 lgeo_.emplace(creator_->localGeometry(element));
103 localValues_.resize(numNodes);
104 for (std::int64_t i = shift, i0 = 0; i < (*offsets_)[insertionIndex]; ++i, ++i0) {
105 std::int64_t idx = (*connectivity_)[i];
106 DynamicVector<Field>& v = localValues_[i0];
108 for (
unsigned int j = 0; j < comp_; ++j)
109 v[j] = (*values_)[comp_*idx + j];
126 auto const& lb = lfe_->localBasis();
127 lb.evaluateFunction(lgeo_->global(local), shapeValues_);
128 assert(shapeValues_.size() == localValues_.size());
130 Range y(comp_, Field(0));
131 for (std::size_t i = 0; i < shapeValues_.size(); ++i)
132 y.axpy(shapeValues_[i], localValues_[i]);
138 GridCreator
const* creator_ =
nullptr;
139 std::vector<Field>
const* values_ =
nullptr;
141 std::vector<std::uint8_t>
const* types_ =
nullptr;
142 std::vector<std::int64_t>
const* offsets_ =
nullptr;
143 std::vector<std::int64_t>
const* connectivity_ =
nullptr;
146 std::optional<LFE> lfe_ = std::nullopt;
147 std::optional<typename GridCreator::LocalGeometry> lgeo_ = std::nullopt;
150 std::vector<DynamicVector<Field>> localValues_;
151 mutable std::vector<typename LB::Traits::RangeType> shapeValues_;
156 class CellDataLocalFunction
159 using LocalContext = LC;
160 using Domain =
typename LC::Geometry::LocalCoordinate;
161 using Range = DynamicVector<Field>;
166 CellDataLocalFunction (GridCreator
const* creator, std::vector<Field>
const* values,
unsigned int comp,
167 std::vector<std::uint8_t>
const* ,
168 std::vector<std::int64_t>
const* ,
169 std::vector<std::int64_t>
const* )
175 CellDataLocalFunction () =
default;
179 void bind (LocalContext
const& element)
181 unsigned int idx = creator_->factory().insertionIndex(element);
184 DynamicVector<Field>& v = localValue_;
187 for (
unsigned int j = 0; j < comp_; ++j)
188 v[j] = (*values_)[comp_*idx + j];
203 GridCreator
const* creator_ =
nullptr;
204 std::vector<Field>
const* values_ =
nullptr;
208 DynamicVector<Field> localValue_;
213 using LocalFunction = std::conditional_t< std::is_same_v<Context,PointContext>,
214 PointDataLocalFunction<LC>,
215 CellDataLocalFunction<LC>>;
222 std::vector<std::uint8_t>
const& types,
223 std::vector<std::int64_t>
const& offsets,
224 std::vector<std::int64_t>
const& connectivity)
227 , name_(std::move(
name))
232 , connectivity_(&connectivity)
240 DUNE_THROW(Dune::NotImplemented,
"Evaluation in global coordinates not implemented.");
241 return Range(ncomps_, 0);
250 std::string
const&
name ()
const
269 return {gf.creator_, gf.values_, gf.ncomps_, gf.types_, gf.offsets_, gf.connectivity_};
273 GridCreator
const* creator_ =
nullptr;
274 std::vector<Field>
const* values_ =
nullptr;
275 std::string name_ =
"GridFunction";
276 unsigned int ncomps_ = 0;
278 std::vector<std::uint8_t>
const* types_ =
nullptr;
279 std::vector<std::int64_t>
const* offsets_ =
nullptr;
280 std::vector<std::int64_t>
const* connectivity_ =
nullptr;
282 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
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
friend LocalFunction< typename EntitySet::Element > localFunction(LagrangeGridFunction const &gf)
Definition: lagrangegridfunction.hh:267
Range(Domain) Signature
Definition: lagrangegridfunction.hh:43
int numComponents() const
Definition: lagrangegridfunction.hh:255
typename EntitySet::GlobalCoordinate Domain
Definition: lagrangegridfunction.hh:41
std::string const & name() const
Definition: lagrangegridfunction.hh:250
Range operator()(Domain const &global) const
Global evaluation. Not supported!
Definition: lagrangegridfunction.hh:238
LagrangeGridFunction()=default
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:220
Vtk::DataTypes dataType() const
Definition: lagrangegridfunction.hh:260
EntitySet const & entitySet() const
Return a type that defines the element that can be iterated.
Definition: lagrangegridfunction.hh:245
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