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_UTILITY_HH 5 #define DUNE_PDELAB_ORDERING_UTILITY_HH 6 7 #include <vector> 8 #include <bitset> 9 10 #include <dune/pdelab/common/dofindex.hh> 11 #include <dune/pdelab/common/globaldofindex.hh> 12 #include <dune/pdelab/ordering/transformations.hh> 13 14 namespace Dune { 15 namespace PDELab { 16 17 //! \addtogroup Ordering 18 //! \{ 19 20 //! Index merging algorithm for global orderings. 21 struct MergeMode 22 { 23 24 enum type { 25 lexicographic, //!< Lexicographically ordered ([i1,i2],[j1,j2] -> [i1,i2,j1,j2]). 26 interleaved //!< Indices are interleaved according to a user-supplied pattern ([i1,i2],[j1,j2] -> [i1,j1,i2,j2]). 27 }; 28 29 }; 30 31 //! Information about order semantics on multi-indices 32 enum class MultiIndexOrder { 33 //! indices are ordered from inner to outer container: {inner,...,outer} 34 Inner2Outer, 35 //! indices are ordered from outer to inner container: {outer,...,inner} 36 Outer2Inner, 37 }; 38 39 40 #ifndef DOXYGEN 41 42 namespace ordering { 43 44 // This is an implementation detail of the composite orderings, no need to confuse our users! 45 struct update_direct_children 46 : public TypeTree::DirectChildrenVisitor 47 , public TypeTree::DynamicTraversal 48 { 49 50 template<typename GFS, typename Child, typename TreePath, typename ChildIndex> afterChildDune::PDELab::ordering::update_direct_children51 void afterChild(const GFS& gfs, Child& child, TreePath tp, ChildIndex childIndex) const 52 { 53 child.update(); 54 } 55 56 }; 57 58 } // end namespace ordering 59 60 #endif // DOXYGEN 61 62 63 struct DefaultDOFIndexAccessor 64 { 65 66 template<typename DOFIndex, typename SizeType, typename IndexType> 67 static typename std::enable_if< 68 std::is_integral<IndexType>::value 69 >::type storeDune::PDELab::DefaultDOFIndexAccessor70 store(DOFIndex& dof_index, const GeometryType& gt, SizeType entity_index, IndexType tree_index) 71 { 72 dof_index.clear(); 73 dof_index.entityIndex()[0] = GlobalGeometryTypeIndex::index(gt); 74 dof_index.entityIndex()[1] = entity_index; 75 dof_index.treeIndex().push_back(tree_index); 76 } 77 78 template<typename DOFIndex, typename SizeType, typename IndexType> 79 static typename std::enable_if< 80 !std::is_integral<IndexType>::value 81 >::type storeDune::PDELab::DefaultDOFIndexAccessor82 store(DOFIndex& dof_index, const GeometryType& gt, SizeType entity_index, IndexType tree_index) 83 { 84 dof_index.entityIndex()[0] = GlobalGeometryTypeIndex::index(gt); 85 dof_index.entityIndex()[1] = entity_index; 86 dof_index.treeIndex() = tree_index; 87 } 88 89 template<typename DOFIndex, typename SizeType, typename IndexType> 90 static typename std::enable_if< 91 std::is_integral<IndexType>::value 92 >::type storeDune::PDELab::DefaultDOFIndexAccessor93 store(DOFIndex& dof_index, SizeType gt_index, SizeType entity_index, IndexType tree_index) 94 { 95 dof_index.clear(); 96 dof_index.entityIndex()[0] = gt_index; 97 dof_index.entityIndex()[1] = entity_index; 98 dof_index.treeIndex().push_back(tree_index); 99 } 100 101 template<typename DOFIndex, typename SizeType, typename IndexType> 102 static typename std::enable_if< 103 !std::is_integral<IndexType>::value 104 >::type storeDune::PDELab::DefaultDOFIndexAccessor105 store(DOFIndex& dof_index, SizeType gt_index, SizeType entity_index, IndexType tree_index) 106 { 107 dof_index.entityIndex()[0] = gt_index; 108 dof_index.entityIndex()[1] = entity_index; 109 dof_index.treeIndex() = tree_index; 110 } 111 112 113 struct GeometryIndex 114 { 115 116 template<typename Index> geometryTypeDune::PDELab::DefaultDOFIndexAccessor::GeometryIndex117 static std::size_t geometryType(const Index& geometry_index) 118 { 119 return geometry_index[0]; 120 } 121 122 template<typename Index> entityIndexDune::PDELab::DefaultDOFIndexAccessor::GeometryIndex123 static std::size_t entityIndex(const Index& geometry_index) 124 { 125 return geometry_index[1]; 126 } 127 128 template<typename Index, typename SizeType> storeDune::PDELab::DefaultDOFIndexAccessor::GeometryIndex129 static void store(Index& index, const GeometryType& gt, SizeType entity_index) 130 { 131 index[0] = GlobalGeometryTypeIndex::index(gt); 132 index[1] = entity_index; 133 } 134 135 template<typename Index, typename SizeType> storeDune::PDELab::DefaultDOFIndexAccessor::GeometryIndex136 static void store(Index& index, SizeType geometry_index, SizeType entity_index) 137 { 138 index[0] = geometry_index; 139 index[1] = entity_index; 140 } 141 142 }; 143 144 template<typename DOFIndex> geometryTypeDune::PDELab::DefaultDOFIndexAccessor145 static std::size_t geometryType(const DOFIndex& dof_index) 146 { 147 return GeometryIndex::geometryType(dof_index.entityIndex()); 148 } 149 150 template<typename DOFIndex> entityIndexDune::PDELab::DefaultDOFIndexAccessor151 static std::size_t entityIndex(const DOFIndex& dof_index) 152 { 153 return GeometryIndex::entityIndex(dof_index.entityIndex()); 154 } 155 156 }; 157 158 struct SimpleDOFIndexAccessor 159 { 160 161 template<typename DOFIndex, typename SizeType> storeDune::PDELab::SimpleDOFIndexAccessor162 static void store(DOFIndex& dof_index, const GeometryType& gt, SizeType entity_index, SizeType tree_index) 163 { 164 dof_index = entity_index; 165 } 166 167 }; 168 169 170 template<typename DI, typename CI, MultiIndexOrder CIOrder = MultiIndexOrder::Inner2Outer> 171 struct SimpleOrderingTraits 172 { 173 174 typedef DI DOFIndex; 175 176 typedef CI ContainerIndex; 177 178 typedef std::size_t SizeType; 179 180 typedef DefaultDOFIndexAccessor DOFIndexAccessor; 181 182 //! Inform about ContainerIndex multi-index order semantics 183 static constexpr MultiIndexOrder ContainerIndexOrder = CIOrder; 184 }; 185 186 187 template<typename SizeType_, typename CI, MultiIndexOrder CIOrder> 188 struct SimpleOrderingTraits<SimpleDOFIndex<SizeType_>,CI,CIOrder> 189 { 190 191 typedef SimpleDOFIndex<SizeType_> DOFIndex; 192 193 typedef CI ContainerIndex; 194 195 typedef SizeType_ SizeType; 196 197 typedef SimpleDOFIndexAccessor DOFIndexAccessor; 198 199 //! Inform about ContainerIndex multi-index order semantics 200 static constexpr MultiIndexOrder ContainerIndexOrder = CIOrder; 201 }; 202 203 template <typename DI, typename CI, 204 MultiIndexOrder CIOrder = MultiIndexOrder::Inner2Outer> 205 struct OrderingTraits : public SimpleOrderingTraits<DI, CI, CIOrder> { 206 207 // The maximum dimension supported (length of bitsets) 208 // 32 dimensions should probably be fine for now... ;-) 209 static const std::size_t max_dim = 32; 210 211 typedef std::bitset<max_dim> CodimFlag; 212 213 typedef typename DI::TreeIndex TreeIndex; 214 215 typedef typename DI::View DOFIndexView; 216 typedef typename DI::View::TreeIndex TreeIndexView; 217 218 typedef typename DI::size_type SizeType; 219 typedef typename DI::size_type size_type; 220 }; 221 222 template <typename ES, typename DI, typename CI, 223 MultiIndexOrder CIOrder = MultiIndexOrder::Outer2Inner> 224 struct LocalOrderingTraits : public OrderingTraits<DI, CI, CIOrder> { 225 226 using EntitySet = ES; 227 using GridView = typename ES::GridView; 228 }; 229 230 template<typename ES, typename DI, typename CI> 231 struct GridViewOrderingTraits 232 : public LocalOrderingTraits<ES,DI,CI> 233 { 234 235 typedef typename DI::EntityIndex EntityIndex; 236 typedef typename DI::View::EntityIndex EntityIndexView; 237 238 }; 239 240 241 template<typename DI, typename CI> 242 class VirtualOrderingBase 243 { 244 public: 245 246 typedef OrderingTraits<DI,CI> Traits; 247 VirtualOrderingBase()248 VirtualOrderingBase() {} 249 virtual ~VirtualOrderingBase() = default; 250 251 virtual void map_index_dynamic(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const = 0; 252 }; 253 254 255 template<typename child_type> 256 struct extract_child_bases 257 : public TypeTree::DirectChildrenVisitor 258 , public TypeTree::DynamicTraversal 259 { 260 261 template<typename Node, typename Child, typename TreePath, typename ChildIndex> afterChildDune::PDELab::extract_child_bases262 void afterChild(const Node& node, Child& child, TreePath tp, ChildIndex child_index) 263 { 264 extract_child(node,child,child_index); 265 } 266 267 template<typename Node, typename Child, typename ChildIndex> 268 typename std::enable_if<Node::has_dynamic_ordering_children>::type extract_childDune::PDELab::extract_child_bases269 extract_child(const Node& node, Child& child, ChildIndex child_index) 270 { 271 _children[child_index] = &child; 272 } 273 274 template<typename Node, typename Child, typename ChildIndex> 275 typename std::enable_if<!Node::has_dynamic_ordering_children>::type extract_childDune::PDELab::extract_child_bases276 extract_child(const Node& node, Child& child, ChildIndex child_index) 277 { 278 } 279 extract_child_basesDune::PDELab::extract_child_bases280 extract_child_bases(std::vector<child_type*>& children) 281 : _children(children) 282 {} 283 284 private: 285 std::vector<child_type*>& _children; 286 287 }; 288 289 290 //! Dummy iterator type over DOF indices. 291 /** 292 * This dummy iterator is used to support omitting the calculation 293 * of DOFIndex values in the per-entity index lookup methods of 294 * orderings. By defining all operations performed on the DOFIndex 295 * iterator and its value by this methods as no-ops, we can reuse the 296 * combined implementation mapping both DOFIndex and ContainerIndex for 297 * the (much more common) case of only having to map the ContainerIndex 298 * values. 299 */ 300 struct DummyDOFIndexIterator 301 { 302 303 typedef std::size_t size_type; 304 operator ++Dune::PDELab::DummyDOFIndexIterator305 DummyDOFIndexIterator& operator++() 306 { 307 return *this; 308 } 309 operator +=Dune::PDELab::DummyDOFIndexIterator310 DummyDOFIndexIterator& operator+=(size_type i) 311 { 312 return *this; 313 } 314 operator *Dune::PDELab::DummyDOFIndexIterator315 DummyDOFIndexIterator& operator*() 316 { 317 return *this; 318 } 319 operator ->Dune::PDELab::DummyDOFIndexIterator320 DummyDOFIndexIterator* operator->() 321 { 322 return this; 323 } 324 treeIndexDune::PDELab::DummyDOFIndexIterator325 DummyDOFIndexIterator& treeIndex() 326 { 327 return *this; 328 } 329 operator ==Dune::PDELab::DummyDOFIndexIterator330 bool operator==(const DummyDOFIndexIterator& r) const 331 { 332 return true; 333 } 334 operator !=Dune::PDELab::DummyDOFIndexIterator335 bool operator!=(const DummyDOFIndexIterator& r) const 336 { 337 return !operator==(r); 338 } 339 push_backDune::PDELab::DummyDOFIndexIterator340 void push_back(size_type i) 341 {} 342 343 }; 344 345 /** 346 * @brief Adapter to create a size provider from an ordering 347 * @details This adapter is meant to be used in allocation and 348 * resizing of vectors containers. 349 * In particular, this adapter is needed because the ordering library give 350 * sizes for multi-indices ordered with Inner2Outer semantis, while 351 * resizing algorithms are faster and easier when using Outer2Inner 352 * semantics. 353 * 354 * - This class makes type erasure on the ordering type. 355 * - This class has value semantics. 356 * - This class always receives container indices with Outer2Inner order 357 * 358 * @tparam SizeType_ return type of the size method 359 * @tparam ContainerIndex_ argument type of the size method 360 * @tparam OriginOrder enum with MultiIndexOrder semantics of the origin ordering 361 */ 362 template <class Size, class ContainerIndex_, MultiIndexOrder OriginOrder> 363 struct SizeProviderAdapter { 364 365 /** 366 * @brief Construct a new Size Provider Adapter object 367 * 368 * @tparam Ordering The type of the ordering to adapt 369 * @param ordering A shared pointer to the ordering 370 */ 371 template <class Ordering> SizeProviderAdapterDune::PDELab::SizeProviderAdapter372 SizeProviderAdapter(const std::shared_ptr<const Ordering> &ordering) 373 : _size_provider([=](const ContainerIndex_ &partial_multiindex) { 374 return ordering->size(partial_multiindex); 375 }) { 376 static_assert(Ordering::Traits::ContainerIndexOrder == OriginOrder); 377 } 378 379 //! Partial MultiIndex of a ContainerIndex 380 using ContainerIndex = ContainerIndex_; 381 382 //! Partial MultiIndex of a ContainerIndex 383 using SizePrefix = ContainerIndex; 384 385 //! Type that refers to the size of containers 386 using SizeType = Size; 387 388 //! Inform about ContainerIndex multi-index order semantics 389 static constexpr MultiIndexOrder ContainerIndexOrder = MultiIndexOrder::Outer2Inner; 390 391 /** 392 * @brief Gives the size for a given prefix 393 * @param prefix MultiIndex with a partial path to a container 394 * @return Traits::SizeType The size required for such a path 395 */ sizeDune::PDELab::SizeProviderAdapter396 SizeType size(const SizePrefix &prefix) const { 397 if constexpr (OriginOrder == MultiIndexOrder::Inner2Outer) { 398 // reversing Outer2Inner prefix into a Inner2Outer suffix 399 ContainerIndex suffix; 400 suffix.resize(prefix.size()); 401 std::reverse_copy(prefix.begin(), prefix.end(),suffix.begin()); 402 // forward size request to ordering with new Inner2Outer suffix 403 return _size_provider(suffix); 404 } else { 405 // prefix is already Outer2Inner, forward to size provider 406 return _size_provider(prefix); 407 } 408 } 409 410 private: 411 412 const std::function<SizeType(const ContainerIndex &)> _size_provider; 413 }; 414 415 //! template deduction guide for orderings 416 template<class Ordering> 417 SizeProviderAdapter(const std::shared_ptr<const Ordering>& ordering) 418 -> SizeProviderAdapter<typename Ordering::Traits::SizeType, 419 typename Ordering::Traits::ContainerIndex, 420 Ordering::Traits::ContainerIndexOrder>; 421 422 423 //! \} group Ordering 424 } // namespace PDELab 425 } // namespace Dune 426 427 #endif // DUNE_PDELAB_ORDERING_UTILITY_HH 428