6#ifndef DUNE_TYPETREE_TREECONTAINER_HH
7#define DUNE_TYPETREE_TREECONTAINER_HH
14#include <dune/common/indices.hh>
15#include <dune/common/hybridutilities.hh>
19#include <dune/typetree/treepath.hh>
37 template<
class LeafToValue>
38 class ContainerFactory
41 using DynamicDegreeConcept =
decltype((std::size_t(std::declval<N>().
degree()),
true));
44 using StaticDegreeConcept =
decltype((std::integral_constant<std::size_t,
N::degree()>{},
true));
47 using DynamicChildAccessConcept =
decltype((std::declval<N>().child(0u),
true));
58 ContainerFactory(LeafToValue leafToValue) :
59 leafToValue_(leafToValue)
63 auto operator()(
const Node& node)
71 std::enable_if_t<Node::isLeaf, bool> =
true>
74 return leafToValue_(node);
78 StaticDegreeConcept<Node> =
true,
79 DynamicChildAccessConcept<Node> =
true>
83 return std::array{(*this)(node.child(indices))...};
84 }, std::make_index_sequence<std::size_t(Node::degree())>());
88 DynamicDegreeConcept<Node> =
true,
89 DynamicChildAccessConcept<Node> =
true>
92 using TransformedChild =
decltype((*this)(node.child(0)));
93 std::vector<TransformedChild> container;
94 container.reserve(node.degree());
95 for (std::size_t i = 0; i < node.degree(); ++i)
96 container.emplace_back((*
this)(node.child(i)));
101 StaticDegreeConcept<Node> =
true>
105 return Dune::makeTupleVector((*
this)(node.child(indices))...);
106 }, std::make_index_sequence<std::size_t(Node::degree())>());
110 LeafToValue leafToValue_;
117 template<
class Container>
118 class TreeContainerVectorBackend
121 static constexpr decltype(
auto) accessByTreePath(C&& container,
const HybridTreePath<>& path)
126 template<
class C,
class... T>
127 static constexpr decltype(
auto) accessByTreePath(C&& container,
const HybridTreePath<T...>& path)
132 }, std::make_index_sequence<
sizeof...(T)-1>());
133 return accessByTreePath(container[
head], tailPath);
136 template<
class C,
class Tree,
137 std::enable_if_t<Tree::isLeaf, bool> =
true>
143 template<
class C,
class Tree,
144 class =
decltype(std::declval<C>().resize(0u))>
147 container.resize(tree.degree());
149 resizeImpl(container[i], tree.child(i), Dune::PriorityTag<5>{});
153 template<
class C,
class Tree>
157 resizeImpl(container[i], tree.child(i), Dune::PriorityTag<5>{});
162 using TypeTreeConcept =
decltype((
163 std::declval<T>().degree(),
171 TreeContainerVectorBackend(Container&& container) :
172 container_(
std::move(container))
176 template <
class Tree, TypeTreeConcept<Tree> = true>
177 TreeContainerVectorBackend(
const Tree& tree) :
178 TreeContainerVectorBackend()
184 template <
class C = Container,
185 std::enable_if_t<std::is_default_constructible_v<C>,
bool> =
true>
186 TreeContainerVectorBackend() :
191 decltype(
auto)
operator[](
const HybridTreePath<T...>& path)
const
193 return accessByTreePath(container_, path);
197 decltype(
auto)
operator[](
const HybridTreePath<T...>& path)
199 return accessByTreePath(container_, path);
203 template<
class Tree, TypeTreeConcept<Tree> = true>
204 void resize(
const Tree& tree)
209 const Container& data()
const
220 Container container_;
223 template<
class Container>
224 auto makeTreeContainerVectorBackend(Container&& container)
226 return TreeContainerVectorBackend<std::decay_t<Container>>(std::forward<Container>(container));
236 template<
template<
class Node>
class LeafToValue>
237 struct LeafToDefaultConstructibleValue
240 auto operator()(
const Node& node)
const
242 return LeafToValue<Node>{};
267 template<
class Tree,
class LeafToValue>
270 auto f = std::ref(leafToValue);
271 auto factory = Detail::ContainerFactory<decltype(f)>(f);
272 return Detail::makeTreeContainerVectorBackend(factory(tree));
290 template<
class Value,
class Tree>
299 template<
class Value,
class Tree>
305 template<
template<
class Node>
class LeafToValue,
class Tree>
306 using TreeContainer = std::decay_t<decltype(makeTreeContainer(std::declval<const Tree&>(), std::declval<Detail::LeafToDefaultConstructibleValue<LeafToValue>>()))>;
constexpr index_constant< 0 > _0
Compile time index with value 0.
Definition: indices.hh:52
decltype(auto) constexpr unpackIntegerSequence(F &&f, std::integer_sequence< I, i... > sequence)
Unpack an std::integer_sequence<I,i...> to std::integral_constant<I,i>...
Definition: indices.hh:124
std::integral_constant< std::size_t, i > index_constant
An index constant with value i.
Definition: indices.hh:29
constexpr void forEach(Range &&range, F &&f)
Range based for loop.
Definition: hybridutilities.hh:256
std::size_t degree(const Node &node)
Returns the degree of node as run time information.
Definition: nodeinterface.hh:79
constexpr auto treePath(const T &... t)
Constructs a new HybridTreePath from the given indices.
Definition: treepath.hh:326
std::decay_t< decltype(makeTreeContainer< Value >(std::declval< const Tree & >()))> UniformTreeContainer
Alias to container type generated by makeTreeContainer for given tree type and uniform value type.
Definition: treecontainer.hh:300
std::decay_t< decltype(makeTreeContainer(std::declval< const Tree & >(), std::declval< Detail::LeafToDefaultConstructibleValue< LeafToValue > >()))> TreeContainer
Alias to container type generated by makeTreeContainer for give tree type and when using LeafToValue ...
Definition: treecontainer.hh:306
auto makeTreeContainer(const Tree &tree)
Create container havin the same structure as the given tree.
Definition: treecontainer.hh:291
Dune namespace.
Definition: alignedallocator.hh:13
constexpr std::integral_constant< T, I0 > head(std::integer_sequence< T, I0, II... >)
For a sequence [head,tail...) return the single head element.
Definition: integersequence.hh:53
Utilities for reduction like operations on ranges.
Helper class for tagging priorities.
Definition: typeutilities.hh:87
Helper class for tagging priorities.
Definition: typeutilities.hh:73
Provides the TupleVector class that augments std::tuple by operator[].