1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- 2 // vi: set et ts=4 sw=2 sts=2: 3 #ifndef DUNE_AMGHIERARCHY_HH 4 #define DUNE_AMGHIERARCHY_HH 5 6 #include <list> 7 #include <memory> 8 #include <limits> 9 #include <dune/common/stdstreams.hh> 10 #include <dune/common/timer.hh> 11 #include <dune/common/bigunsignedint.hh> 12 #include <dune/istl/paamg/construction.hh> 13 14 namespace Dune 15 { 16 namespace Amg 17 { 18 /** 19 * @addtogroup ISTL_PAAMG 20 * 21 * @{ 22 */ 23 24 /** @file 25 * @author Markus Blatt 26 * @brief Provides a classes representing the hierarchies in AMG. 27 */ 28 29 /** 30 * @brief A hierarchy of containers (e.g. matrices or vectors) 31 * 32 * Because sometimes a redistribution of the parallel data might be 33 * advisable one can add redistributed version of the container at 34 * each level. 35 */ 36 template<typename T, typename A=std::allocator<T> > 37 class Hierarchy 38 { 39 public: 40 /** 41 * @brief The type of the container we store. 42 */ 43 typedef T MemberType; 44 45 template<typename T1, typename T2> 46 class LevelIterator; 47 48 private: 49 /** 50 * @brief An element in the hierarchy. 51 */ 52 struct Element 53 { 54 friend class LevelIterator<Hierarchy<T,A>, T>; 55 friend class LevelIterator<const Hierarchy<T,A>, const T>; 56 57 /** @brief The next coarser element in the list. */ 58 std::weak_ptr<Element> coarser_; 59 60 /** @brief The next finer element in the list. */ 61 std::shared_ptr<Element> finer_; 62 63 /** @brief Pointer to the element. */ 64 std::shared_ptr<MemberType> element_; 65 66 /** @brief The redistributed version of the element. */ 67 std::shared_ptr<MemberType> redistributed_; 68 }; 69 public: 70 71 /** 72 * @brief The allocator to use for the list elements. 73 */ 74 using Allocator = typename std::allocator_traits<A>::template rebind_alloc<Element>; 75 76 typedef typename ConstructionTraits<T>::Arguments Arguments; 77 78 /** 79 * @brief Construct a new hierarchy. 80 * @param first std::shared_ptr to the first element in the hierarchy. 81 */ 82 Hierarchy(const std::shared_ptr<MemberType> & first); 83 84 /** 85 * @brief Construct an empty hierarchy. 86 */ Hierarchy()87 Hierarchy() : levels_(0) 88 {} 89 90 /** 91 * @brief Copy constructor (deep copy!). 92 */ 93 Hierarchy(const Hierarchy& other); 94 95 /** 96 * @brief Add an element on a coarser level. 97 * @param args The arguments needed for the construction. 98 */ 99 void addCoarser(Arguments& args); 100 101 void addRedistributedOnCoarsest(Arguments& args); 102 103 /** 104 * @brief Add an element on a finer level. 105 * @param args The arguments needed for the construction. 106 */ 107 void addFiner(Arguments& args); 108 109 /** 110 * @brief Iterator over the levels in the hierarchy. 111 * 112 * operator++() moves to the next coarser level in the hierarchy. 113 * while operator--() moves to the next finer level in the hierarchy. 114 */ 115 template<class C, class T1> 116 class LevelIterator 117 : public BidirectionalIteratorFacade<LevelIterator<C,T1>,T1,T1&> 118 { 119 friend class LevelIterator<typename std::remove_const<C>::type, 120 typename std::remove_const<T1>::type >; 121 friend class LevelIterator<const typename std::remove_const<C>::type, 122 const typename std::remove_const<T1>::type >; 123 124 public: 125 /** @brief Constructor. */ LevelIterator()126 LevelIterator() 127 {} 128 LevelIterator(std::shared_ptr<Element> element)129 LevelIterator(std::shared_ptr<Element> element) 130 : element_(element) 131 {} 132 133 /** @brief Copy constructor. */ LevelIterator(const LevelIterator<typename std::remove_const<C>::type,typename std::remove_const<T1>::type> & other)134 LevelIterator(const LevelIterator<typename std::remove_const<C>::type, 135 typename std::remove_const<T1>::type>& other) 136 : element_(other.element_) 137 {} 138 139 /** @brief Copy constructor. */ LevelIterator(const LevelIterator<const typename std::remove_const<C>::type,const typename std::remove_const<T1>::type> & other)140 LevelIterator(const LevelIterator<const typename std::remove_const<C>::type, 141 const typename std::remove_const<T1>::type>& other) 142 : element_(other.element_) 143 {} 144 145 /** 146 * @brief Equality check. 147 */ equals(const LevelIterator<typename std::remove_const<C>::type,typename std::remove_const<T1>::type> & other) const148 bool equals(const LevelIterator<typename std::remove_const<C>::type, 149 typename std::remove_const<T1>::type>& other) const 150 { 151 return element_ == other.element_; 152 } 153 154 /** 155 * @brief Equality check. 156 */ equals(const LevelIterator<const typename std::remove_const<C>::type,const typename std::remove_const<T1>::type> & other) const157 bool equals(const LevelIterator<const typename std::remove_const<C>::type, 158 const typename std::remove_const<T1>::type>& other) const 159 { 160 return element_ == other.element_; 161 } 162 163 /** @brief Dereference the iterator. */ dereference() const164 T1& dereference() const 165 { 166 return *(element_->element_); 167 } 168 169 /** @brief Move to the next coarser level */ increment()170 void increment() 171 { 172 element_ = element_->coarser_.lock(); 173 } 174 175 /** @brief Move to the next fine level */ decrement()176 void decrement() 177 { 178 element_ = element_->finer_; 179 } 180 181 /** 182 * @brief Check whether there was a redistribution at the current level. 183 * @return True if there is a redistributed version of the container at the current level. 184 */ isRedistributed() const185 bool isRedistributed() const 186 { 187 return (bool)element_->redistributed_; 188 } 189 190 /** 191 * @brief Get the redistributed container. 192 * @return The redistributed container. 193 */ getRedistributed() const194 T1& getRedistributed() const 195 { 196 assert(element_->redistributed_); 197 return *element_->redistributed_; 198 } addRedistributed(std::shared_ptr<T1> t)199 void addRedistributed(std::shared_ptr<T1> t) 200 { 201 element_->redistributed_ = t; 202 } 203 deleteRedistributed()204 void deleteRedistributed() 205 { 206 element_->redistributed_ = nullptr; 207 } 208 209 private: 210 std::shared_ptr<Element> element_; 211 }; 212 213 /** @brief Type of the mutable iterator. */ 214 typedef LevelIterator<Hierarchy<T,A>,T> Iterator; 215 216 /** @brief Type of the const iterator. */ 217 typedef LevelIterator<const Hierarchy<T,A>, const T> ConstIterator; 218 219 /** 220 * @brief Get an iterator positioned at the finest level. 221 * @return An iterator positioned at the finest level. 222 */ 223 Iterator finest(); 224 225 /** 226 * @brief Get an iterator positioned at the coarsest level. 227 * @return An iterator positioned at the coarsest level. 228 */ 229 Iterator coarsest(); 230 231 232 /** 233 * @brief Get an iterator positioned at the finest level. 234 * @return An iterator positioned at the finest level. 235 */ 236 ConstIterator finest() const; 237 238 /** 239 * @brief Get an iterator positioned at the coarsest level. 240 * @return An iterator positioned at the coarsest level. 241 */ 242 ConstIterator coarsest() const; 243 244 /** 245 * @brief Get the number of levels in the hierarchy. 246 * @return The number of levels. 247 */ 248 std::size_t levels() const; 249 250 private: 251 /** @brief fix memory management of the finest element in the hierarchy 252 253 This object is passed in from outside, so that we have to 254 deal with shared ownership. 255 */ 256 std::shared_ptr<MemberType> originalFinest_; 257 /** @brief The finest element in the hierarchy. */ 258 std::shared_ptr<Element> finest_; 259 /** @brief The coarsest element in the hierarchy. */ 260 std::shared_ptr<Element> coarsest_; 261 /** @brief The allocator for the list elements. */ 262 Allocator allocator_; 263 /** @brief The number of levels in the hierarchy. */ 264 int levels_; 265 }; 266 267 template<class T, class A> Hierarchy(const std::shared_ptr<MemberType> & first)268 Hierarchy<T,A>::Hierarchy(const std::shared_ptr<MemberType> & first) 269 : originalFinest_(first) 270 { 271 finest_ = std::allocate_shared<Element>(allocator_); 272 finest_->element_ = originalFinest_; 273 coarsest_ = finest_; 274 levels_ = 1; 275 } 276 277 //! \brief deep copy of a given hierarchy 278 //TODO: do we actually want to support this? This might be very expensive?! 279 template<class T, class A> Hierarchy(const Hierarchy & other)280 Hierarchy<T,A>::Hierarchy(const Hierarchy& other) 281 : allocator_(other.allocator_), 282 levels_(other.levels_) 283 { 284 if(!other.finest_) 285 { 286 finest_=coarsest_=nullptr; 287 return; 288 } 289 finest_ = std::allocate_shared<Element>(allocator_); 290 std::shared_ptr<Element> finer_; 291 std::shared_ptr<Element> current_ = finest_; 292 std::weak_ptr<Element> otherWeak_ = other.finest_; 293 294 while(! otherWeak_.expired()) 295 { 296 // create shared_ptr from weak_ptr, we just checked that this is safe 297 std::shared_ptr<Element> otherCurrent_ = std::shared_ptr<Element>(otherWeak_); 298 // clone current level 299 //TODO: should we use the allocator? 300 current_->element_ = 301 std::make_shared<MemberType>(*(otherCurrent_->element_)); 302 current_->finer_=finer_; 303 if(otherCurrent_->redistributed_) 304 current_->redistributed_ = 305 std::make_shared<MemberType>(*(otherCurrent_->redistributed_)); 306 finer_=current_; 307 if(not otherCurrent_->coarser_.expired()) 308 { 309 auto c = std::allocate_shared<Element>(allocator_); 310 current_->coarser_ = c; 311 current_ = c; 312 } 313 // go to coarser level 314 otherWeak_ = otherCurrent_->coarser_; 315 } 316 coarsest_=current_; 317 } 318 319 template<class T, class A> levels() const320 std::size_t Hierarchy<T,A>::levels() const 321 { 322 return levels_; 323 } 324 325 template<class T, class A> addRedistributedOnCoarsest(Arguments & args)326 void Hierarchy<T,A>::addRedistributedOnCoarsest(Arguments& args) 327 { 328 coarsest_->redistributed_ = ConstructionTraits<MemberType>::construct(args); 329 } 330 331 template<class T, class A> addCoarser(Arguments & args)332 void Hierarchy<T,A>::addCoarser(Arguments& args) 333 { 334 if(!coarsest_) { 335 // we have no levels at all... 336 assert(!finest_); 337 // allocate into the shared_ptr 338 originalFinest_ = ConstructionTraits<MemberType>::construct(args); 339 coarsest_ = std::allocate_shared<Element>(allocator_); 340 coarsest_->element_ = originalFinest_; 341 finest_ = coarsest_; 342 }else{ 343 auto old_coarsest = coarsest_; 344 coarsest_ = std::allocate_shared<Element>(allocator_); 345 coarsest_->finer_ = old_coarsest; 346 coarsest_->element_ = ConstructionTraits<MemberType>::construct(args); 347 old_coarsest->coarser_ = coarsest_; 348 } 349 ++levels_; 350 } 351 352 353 template<class T, class A> addFiner(Arguments & args)354 void Hierarchy<T,A>::addFiner(Arguments& args) 355 { 356 //TODO: wouldn't it be better to do this in the constructor?' 357 if(!finest_) { 358 // we have no levels at all... 359 assert(!coarsest_); 360 // allocate into the shared_ptr 361 originalFinest_ = ConstructionTraits<MemberType>::construct(args); 362 finest_ = std::allocate_shared<Element>(allocator_); 363 finest_->element = originalFinest_; 364 coarsest_ = finest_; 365 }else{ 366 finest_->finer_ = std::allocate_shared<Element>(allocator_); 367 finest_->finer_->coarser_ = finest_; 368 finest_ = finest_->finer_; 369 finest_->element = ConstructionTraits<T>::construct(args); 370 } 371 ++levels_; 372 } 373 374 template<class T, class A> finest()375 typename Hierarchy<T,A>::Iterator Hierarchy<T,A>::finest() 376 { 377 return Iterator(finest_); 378 } 379 380 template<class T, class A> coarsest()381 typename Hierarchy<T,A>::Iterator Hierarchy<T,A>::coarsest() 382 { 383 return Iterator(coarsest_); 384 } 385 386 template<class T, class A> finest() const387 typename Hierarchy<T,A>::ConstIterator Hierarchy<T,A>::finest() const 388 { 389 return ConstIterator(finest_); 390 } 391 392 template<class T, class A> coarsest() const393 typename Hierarchy<T,A>::ConstIterator Hierarchy<T,A>::coarsest() const 394 { 395 return ConstIterator(coarsest_); 396 } 397 /** @} */ 398 } // namespace Amg 399 } // namespace Dune 400 401 #endif 402