dune-pdelab  2.5-dev
subordering.hh
Go to the documentation of this file.
1 // -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=8 sw=2 sts=2:
3 
4 #ifndef DUNE_PDELAB_ORDERING_SUBORDERING_HH
5 #define DUNE_PDELAB_ORDERING_SUBORDERING_HH
6 
7 #include <algorithm>
8 #include <iterator>
9 #include <memory>
10 
11 #include <dune/common/array.hh>
12 
13 #include <dune/typetree/treepath.hh>
14 #include <dune/typetree/proxynode.hh>
15 #include <dune/typetree/childextraction.hh>
16 
18 
19 namespace Dune {
20  namespace PDELab {
21 
24 
26 
48  template<typename BaseOrdering_, typename TreePath>
50  : public TypeTree::ProxyNode<const TypeTree::ChildForTreePath<BaseOrdering_,TreePath>>
51  {
52 
53  using NodeT = TypeTree::ProxyNode<const TypeTree::ChildForTreePath<BaseOrdering_,TreePath>>;
54 
55  public:
56 
58  typedef BaseOrdering_ BaseOrdering;
59 
61  using TargetOrdering = TypeTree::ChildForTreePath<BaseOrdering,TreePath>;
62 
64  typedef typename BaseOrdering::Traits Traits;
65 
67  typedef typename BaseOrdering::ContainerAllocationTag ContainerAllocationTag;
68 
70  typedef typename BaseOrdering::CacheTag CacheTag;
71 
73  static const bool has_dynamic_ordering_children = TargetOrdering::has_dynamic_ordering_children;
74 
76  static const bool consume_tree_index = TargetOrdering::consume_tree_index;
77 
78 
80 
93  explicit SubOrdering(std::shared_ptr<const BaseOrdering> base_ordering)
94  : NodeT(base_ordering->child(TreePath()))
95  , _base_ordering(base_ordering)
96  {
97  update();
98  }
99 
101  void update()
102  {}
103 
104 
105  template<typename ItIn, typename ItOut>
106  void map_lfs_indices(ItIn begin, const ItIn end, ItOut out) const
107  {
108  // Do the mapping up to the root ordering.
109  // Avoid spelling out the type of ItIn here (it has to be a DOFIndexViewIterator),
110  // so we don't need to include lfsindexcache.hh.
111  map_lfs_indices_to_root_space(TreePath(),
112  begin,
113  end,
114  out);
115  }
116 
117 
118  private:
119 
120 
122  template<typename TP, typename ItIn, typename ItOut>
123  void map_lfs_indices_in_ancestor(TP tp, ItIn& begin, ItIn& end, ItOut out) const
124  {
125  using Ordering = TypeTree::ChildForTreePath<BaseOrdering,TP>;
126 
127  // This logic needs to be replicated from the IndexCache visitor, as we bypass
128  // the tree-visiting algorithm and work our way up the tree all by ourselves.
129  // Don't consume the first entry in the tree path to the parent before it has
130  // been used!
131  if (!std::is_same<TreePath,TP>::value && Ordering::consume_tree_index)
132  {
133  begin.restore_back();
134  end.restore_back();
135  }
136 
137  // Call the single-level mapping step of our ancestor ordering.
138  TypeTree::child(baseOrdering(),tp).map_lfs_indices(begin,end,out);
139  }
140 
141  // Template recursion for walking up the TreePath to the BaseOrdering
142  template<typename TP, typename ItIn, typename ItOut>
143  typename std::enable_if<
145  >::type
146  map_lfs_indices_to_root_space(TP, ItIn begin, ItIn end, ItOut out) const
147  {
148  map_lfs_indices_in_ancestor(TP(),begin,end,out);
149  // recurse further up to the tree
150  map_lfs_indices_to_root_space(typename TypeTree::TreePathPopBack<TP>::type(),begin,end,out);
151  }
152 
153  // End of template recursion for walking up the TreePath to the BaseOrdering
154  template<typename TP, typename ItIn, typename ItOut>
155  typename std::enable_if<
157  >::type
158  map_lfs_indices_to_root_space(TP, ItIn begin, ItIn end, ItOut out) const
159  {
160  map_lfs_indices_in_ancestor(TP(),begin,end,out);
161  }
162 
163  public:
164 
166  typename Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex& di) const
167  {
168  typename Traits::ContainerIndex ci;
169  mapIndex(di,ci);
170  return ci;
171  }
172 
174  void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
175  {
176  baseOrdering().mapIndex(di,ci);
177  }
178 
180  typename Traits::SizeType size() const
181  {
182  return baseOrdering().size();
183  }
184 
186  typename Traits::SizeType blockCount() const
187  {
188  return baseOrdering().blockCount();
189  }
190 
192  typename Traits::SizeType maxLocalSize() const
193  {
194  return targetOrdering().maxLocalSize();
195  }
196 
198  bool contains(typename Traits::SizeType codim) const
199  {
200  return targetOrdering().contains(codim);
201  }
202 
204  bool fixedSize(typename Traits::SizeType codim) const
205  {
206  return targetOrdering().fixedSize(codim);
207  }
208 
210  const BaseOrdering& baseOrdering() const
211  {
212  return *_base_ordering;
213  }
214 
217  {
218  return this->proxiedNode();
219  }
220 
221  private:
222 
223  std::shared_ptr<const BaseOrdering> _base_ordering;
224 
225  };
226 
228 
229  } // namespace PDELab
230 } // namespace Dune
231 
232 #endif // DUNE_PDELAB_ORDERING_SUBORDERING_HH
BaseOrdering::ContainerAllocationTag ContainerAllocationTag
Forwarded tag from BaseOrdering, required by PDELab internals.
Definition: subordering.hh:67
BaseOrdering::CacheTag CacheTag
Forwarded tag from BaseOrdering, required by PDELab internals.
Definition: subordering.hh:70
const TargetOrdering & targetOrdering() const
Returns the TargetOrdering.
Definition: subordering.hh:216
BaseOrdering_ BaseOrdering
The type of the BaseOrdering for which to represent a SubOrdering view.
Definition: subordering.hh:58
BaseOrdering::Traits Traits
Forwarded Ordering traits from BaseOrdering.
Definition: subordering.hh:64
Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex &di) const
Maps di from the DOFIndex subtree to the ContainerIndex in the BaseOrdering.
Definition: subordering.hh:166
Traits::SizeType blockCount() const
Returns the block count of the BaseOrdering.
Definition: subordering.hh:186
void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex &ci) const
Maps di from the DOFIndex subtree to the ContainerIndex in the BaseOrdering - inplace version...
Definition: subordering.hh:174
SubOrdering(std::shared_ptr< const BaseOrdering > base_ordering)
Constructs a SubOrdering for base_ordering.
Definition: subordering.hh:93
TypeTree::ChildForTreePath< BaseOrdering, typename find_ordering_treepath_for_sub_gfs< typename GFS::Ordering, GFS, TreePath >::type > TargetOrdering
The target ordering that makes up the root of this SubOrdering view.
Definition: subordering.hh:61
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
Traits::SizeType size() const
Returns the size of the BaseOrdering.
Definition: subordering.hh:180
A view on a subtree of a multi-component ordering.
Definition: subordering.hh:49
static const bool consume_tree_index
Forwarded ordering property from TargetOrdering, required by PDELab internals.
Definition: subordering.hh:76
const BaseOrdering & baseOrdering() const
Returns the BaseOrdering.
Definition: subordering.hh:210
bool fixedSize(typename Traits::SizeType codim) const
Returns whether the TargetOrdering is of fixed size for entities of codimension codim.
Definition: subordering.hh:204
void map_lfs_indices(ItIn begin, const ItIn end, ItOut out) const
Definition: subordering.hh:106
static const bool has_dynamic_ordering_children
Forwarded ordering property from TargetOrdering, required by PDELab internals.
Definition: subordering.hh:73
void update()
Updates this SubOrdering.
Definition: subordering.hh:101
Traits::SizeType maxLocalSize() const
Returns the maximum per-entity size of the TargetOrdering.
Definition: subordering.hh:192
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
bool contains(typename Traits::SizeType codim) const
Returns whether the TargetOrdering has DOFs attached to entities of codimension codim.
Definition: subordering.hh:198