9 #include <dune/common/exceptions.hh>
10 #include <dune/common/hybridutilities.hh>
11 #include <dune/geometry/utility/typefromvertexcount.hh>
12 #include <dune/geometry/multilineargeometry.hh>
13 #include <dune/localfunctions/lagrange.hh>
14 #include <dune/grid/common/gridfactory.hh>
37 template <
class Gr
idType>
45 using Nodes = std::vector<GlobalCoordinate>;
50 std::vector<std::int64_t>
nodes;
55 using Element =
typename GridType::template Codim<0>::Entity;
62 using LocalGeometry = MultiLinearGeometry<typename Element::Geometry::ctype,Element::dimension,Element::dimension>;
70 std::vector<std::uint64_t>
const& )
82 std::vector<std::int64_t>
const& offsets,
83 std::vector<std::int64_t>
const& connectivity)
85 assert(nodes_.size() > 0);
88 std::vector<std::int64_t> elementVertices(nodes_.size(), -1);
89 parametrization_.reserve(types.size());
91 std::int64_t vertexIndex = 0;
92 for (std::size_t i = 0; i < types.size(); ++i) {
94 if (type.dim() != GridType::dimension)
98 auto refElem = referenceElement<double,GridType::dimension>(type);
100 std::int64_t shift = (i == 0 ? 0 : offsets[i-1]);
101 int nNodes = offsets[i] - shift;
102 int nVertices = refElem.size(GridType::dimension);
105 std::vector<unsigned int> element(nVertices);
106 for (
int j = 0; j < nVertices; ++j) {
107 auto index = connectivity.at(shift + j);
108 auto& vertex = elementVertices.at(index);
110 factory().insertVertex(nodes_.at(index));
111 vertex = vertexIndex++;
117 if (!cellType.noPermutation()) {
119 std::vector<unsigned int> cell(element.size());
120 for (std::size_t j = 0; j < element.size(); ++j)
121 cell[j] = element[cellType.permutation(j)];
122 std::swap(element, cell);
127 auto& param = parametrization_.back();
129 param.nodes.resize(nNodes);
130 for (
int j = 0; j < nNodes; ++j)
131 param.nodes[j] = connectivity.at(shift + j);
132 param.corners = element;
137 factory().insertElement(type, element,
139 }
catch (Dune::GridError
const& ) {
140 factory().insertElement(type, element);
143 factory().insertElement(type, element);
157 assert(!nodes_.empty() && !parametrization_.empty());
158 auto const& localParam = parametrization_.at(insertionIndex);
170 VTK_ASSERT(!nodes_.empty() && !parametrization_.empty());
172 unsigned int insertionIndex =
factory().insertionIndex(element);
173 auto const& localParam = parametrization_.at(insertionIndex);
174 VTK_ASSERT(element.type() == localParam.type);
189 VTK_ASSERT(!nodes_.empty() && !parametrization_.empty());
191 unsigned int insertionIndex =
factory().insertionIndex(element);
192 auto const& localParam = parametrization_.at(insertionIndex);
193 VTK_ASSERT(element.type() == localParam.type);
203 std::vector<unsigned int> indices(element.subEntities(GridType::dimension));
204 for (
unsigned int i = 0; i < element.subEntities(GridType::dimension); ++i)
205 indices[i] =
factory().insertionIndex(element.template subEntity<GridType::dimension>(i));
208 std::vector<unsigned int> permutation(indices.size());
209 for (std::size_t i = 0; i < indices.size(); ++i) {
210 auto it = std::find(localParam.corners.begin(), localParam.corners.end(), indices[i]);
212 permutation[i] = std::distance(localParam.corners.begin(), it);
215 auto refElem = referenceElement<typename Element::Geometry::ctype,Element::dimension>(localParam.type);
216 std::vector<LocalCoordinate> corners(permutation.size());
217 for (std::size_t i = 0; i < permutation.size(); ++i)
218 corners[i] = refElem.position(permutation[i], Element::dimension);
220 return {localParam.type, corners};
226 int order (GeometryType type, std::size_t nNodes)
const
228 for (
int o = 1; o <= int(nNodes); ++o)
229 #
if DUNE_VERSION_LT(DUNE_LOCALFUNCTIONS,2,8)
230 if (numLagrangePoints(type.id(), type.dim(), o) == std::size_t(nNodes))
232 if (numLagrangePoints(type, o) == std::size_t(nNodes))
247 assert(!parametrization_.empty());
248 auto const& localParam = parametrization_.front();
249 return order(localParam);
279 DUNE_THROW(Dune::Exception,
"Cannot pass temporary LagrangeGridCreator to localFunction(). Pass an lvalue-reference instead.");
285 return "LagrangeParametrization";
290 return GlobalCoordinate::size();
295 return dataTypeOf<typename Element::Geometry::ctype>();
307 assert(
false &&
"Should not be used!");
314 assert(
false &&
"Should not be used!");
327 template <
class Gr
id>
331 template <
class Gr
idType,
class FieldType,
class Context>
337 template <
class Gr
id>
340 using ctype =
typename Grid::ctype;
342 using GlobalCoordinate =
typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
343 using LocalCoordinate =
typename Grid::template Codim<0>::Entity::Geometry::LocalCoordinate;
344 using LocalGeometry = MultiLinearGeometry<ctype,Grid::dimension,Grid::dimension>;
346 using LocalFE = LagrangeLocalFiniteElement<Vtk::LagrangePointSet, Grid::dimension, ctype, ctype>;
347 using LocalBasis =
typename LocalFE::Traits::LocalBasisType;
348 using LocalBasisTraits =
typename LocalBasis::Traits;
352 template <
class Nodes,
class LocalParam>
354 : localFE_(param.type,
order)
355 , localNodes_(param.nodes.size())
357 for (std::size_t i = 0; i < localNodes_.size(); ++i)
358 localNodes_[i] = nodes[param.nodes[i]];
362 template <
class Nodes,
class LocalParam,
class LG>
370 template <
class LocalCoordinate>
371 GlobalCoordinate
operator() (LocalCoordinate
const& local)
const
374 LocalCoordinate x = localGeometry_ ? localGeometry_->global(local) : local;
376 LocalBasis
const& localBasis = localFE_.localBasis();
377 localBasis.evaluateFunction(x, shapeValues_);
378 assert(shapeValues_.size() == localNodes_.size());
380 using field_type =
typename LocalBasisTraits::RangeType::field_type;
382 GlobalCoordinate out(0);
383 for (std::size_t i = 0; i < shapeValues_.size(); ++i)
384 out.axpy(field_type(shapeValues_[i]), localNodes_[i]);
391 std::vector<GlobalCoordinate> localNodes_;
392 std::optional<LocalGeometry> localGeometry_;
394 mutable std::vector<typename LocalBasisTraits::RangeType> shapeValues_;
398 template <
class Gr
id>
401 using ctype =
typename Grid::ctype;
402 using LocalContext =
typename Grid::template Codim<0>::Entity;
403 using GlobalCoordinate =
typename LocalContext::Geometry::GlobalCoordinate;
404 using LocalCoordinate =
typename LocalContext::Geometry::LocalCoordinate;
409 : gridCreator_(&gridCreator)
415 void bind (LocalContext
const& element)
417 localContext_ = element;
418 localParametrization_.emplace(gridCreator_->localParametrization(element));
424 GlobalCoordinate
operator() (LocalCoordinate
const& local)
const
426 assert(!!localParametrization_);
427 return (*localParametrization_)(local);
433 return localContext_;
439 LocalContext localContext_;
440 std::optional<LocalParametrization> localParametrization_;
#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
GeometryType to_geometry(std::uint8_t cell)
Definition: types.cc:146
Base class for grid creators in a CRTP style.
Definition: gridcreatorinterface.hh:25
typename Grid::template Codim< 0 >::Entity::Geometry::GlobalCoordinate GlobalCoordinate
Definition: gridcreatorinterface.hh:28
GridType Grid
Definition: gridcreatorinterface.hh:27
GridFactory< Grid > & factory()
Return the associated GridFactory.
Definition: gridcreatorinterface.hh:77
Definition: lagrangegridcreator.hh:40
typename GridType::template Codim< 0 >::Entity Element
Definition: lagrangegridcreator.hh:55
void insertElementsImpl(std::vector< std::uint8_t > const &types, std::vector< std::int64_t > const &offsets, std::vector< std::int64_t > const &connectivity)
Implementation of the interface function insertElements()
Definition: lagrangegridcreator.hh:81
MultiLinearGeometry< typename Element::Geometry::ctype, Element::dimension, Element::dimension > LocalGeometry
Definition: lagrangegridcreator.hh:62
int numComponents() const
Definition: lagrangegridcreator.hh:288
EntitySet entitySet() const
Dummy function returning a placeholder entityset.
Definition: lagrangegridcreator.hh:305
void insertVerticesImpl(std::vector< GlobalCoordinate > const &points, std::vector< std::uint64_t > const &)
Implementation of the interface function insertVertices()
Definition: lagrangegridcreator.hh:69
decltype(std::declval< F >().insertElement(std::declval< GeometryType >(), std::declval< std::vector< unsigned int > const & >(), std::declval< std::function< GlobalCoordinate(LocalCoordinate)> >())) HasParametrizedElements
Definition: lagrangegridcreator.hh:78
friend LocalFunction localFunction(LagrangeGridCreator &&gridCreator)
Definition: lagrangegridcreator.hh:277
int order(GeometryType type, std::size_t nNodes) const
Determine lagrange order from number of points.
Definition: lagrangegridcreator.hh:226
int order() const
Determine lagrange order from number of points from the first element parametrization.
Definition: lagrangegridcreator.hh:245
LocalGeometry localGeometry(Element const &element) const
Construct a transformation of local element coordinates.
Definition: lagrangegridcreator.hh:187
GlobalCoordinate operator()(GlobalCoordinate const &) const
Dummy function returning a placeholder entityset.
Definition: lagrangegridcreator.hh:312
std::string name() const
Definition: lagrangegridcreator.hh:283
int order(ElementParametrization const &localParam) const
Definition: lagrangegridcreator.hh:239
typename Element::Geometry::LocalCoordinate LocalCoordinate
Definition: lagrangegridcreator.hh:56
LocalParametrization localParametrization(Element const &element) const
Construct an element parametrization.
Definition: lagrangegridcreator.hh:168
LocalParametrization localParametrization(unsigned int insertionIndex) const
Construct an element parametrization.
Definition: lagrangegridcreator.hh:155
std::vector< GlobalCoordinate > Nodes
Definition: lagrangegridcreator.hh:45
friend LocalFunction localFunction(LagrangeGridCreator &gridCreator)
Local function representing the parametrization of the grid.
Definition: lagrangegridcreator.hh:267
GridFactory< Grid > & factory()
Return the associated GridFactory.
Definition: gridcreatorinterface.hh:77
Vtk::DataTypes dataType() const
Definition: lagrangegridcreator.hh:293
std::vector< ElementParametrization > Parametrization
Definition: lagrangegridcreator.hh:54
typename Super::GlobalCoordinate GlobalCoordinate
Definition: lagrangegridcreator.hh:43
friend LocalFunction localFunction(LagrangeGridCreator const &gridCreator)
Definition: lagrangegridcreator.hh:272
Definition: lagrangegridcreator.hh:48
std::vector< unsigned int > corners
Definition: lagrangegridcreator.hh:51
GeometryType type
Definition: lagrangegridcreator.hh:49
std::vector< std::int64_t > nodes
Definition: lagrangegridcreator.hh:50
Definition: lagrangegridcreator.hh:299
typename Self::GlobalCoordinate GlobalCoordinate
Definition: lagrangegridcreator.hh:301
GridType Grid
Definition: lagrangegridcreator.hh:300
Definition: lagrangegridcreator.hh:339
LocalParametrization(Nodes const &nodes, LocalParam const ¶m, int order)
Construct a local element parametrization.
Definition: lagrangegridcreator.hh:353
LocalParametrization(Nodes const &nodes, LocalParam const ¶m, int order, LG &&localGeometry)
Construct a local element parametrization for elements with permuted corners.
Definition: lagrangegridcreator.hh:363
Definition: lagrangegridcreator.hh:400
LocalFunction(LagrangeGridCreator &&gridCreator)=delete
LocalFunction(LagrangeGridCreator const &gridCreator)
Definition: lagrangegridcreator.hh:408
void unbind()
Definition: lagrangegridcreator.hh:421
void bind(LocalContext const &element)
Collect a local parametrization on the element.
Definition: lagrangegridcreator.hh:415
LocalContext const & localContext() const
Return the bound element.
Definition: lagrangegridcreator.hh:431
Type-Traits to associate a GridFunction to a GridCreator.
Definition: gridfunctions/common.hh:26
Grid-function representing values from a VTK file with local Lagrange interpolation of the values sto...
Definition: lagrangegridfunction.hh:25
Mapping of Dune geometry types to VTK cell types.
Definition: types.hh:160