1 // 2 // Copyright (c) 2000-2002 3 // Joerg Walter, Mathias Koch 4 // 5 // Distributed under the Boost Software License, Version 1.0. (See 6 // accompanying file LICENSE_1_0.txt or copy at 7 // http://www.boost.org/LICENSE_1_0.txt) 8 // 9 // The authors gratefully acknowledge the support of 10 // GeNeSys mbH & Co. KG in producing this work. 11 // 12 13 #ifndef _BOOST_UBLAS_VECTOR_SPARSE_ 14 #define _BOOST_UBLAS_VECTOR_SPARSE_ 15 16 #include <boost/numeric/ublas/storage_sparse.hpp> 17 #include <boost/numeric/ublas/vector_expression.hpp> 18 #include <boost/numeric/ublas/detail/vector_assign.hpp> 19 #if BOOST_UBLAS_TYPE_CHECK 20 #include <boost/numeric/ublas/vector.hpp> 21 #endif 22 23 // Iterators based on ideas of Jeremy Siek 24 25 namespace boost { namespace numeric { namespace ublas { 26 27 #ifdef BOOST_UBLAS_STRICT_VECTOR_SPARSE 28 29 template<class V> 30 class sparse_vector_element: 31 public container_reference<V> { 32 public: 33 typedef V vector_type; 34 typedef typename V::size_type size_type; 35 typedef typename V::value_type value_type; 36 typedef const value_type &const_reference; 37 typedef value_type *pointer; 38 39 private: 40 // Proxied element operations get_d() const41 void get_d () const { 42 pointer p = (*this) ().find_element (i_); 43 if (p) 44 d_ = *p; 45 else 46 d_ = value_type/*zero*/(); 47 } 48 set(const value_type & s) const49 void set (const value_type &s) const { 50 pointer p = (*this) ().find_element (i_); 51 if (!p) 52 (*this) ().insert_element (i_, s); 53 else 54 *p = s; 55 } 56 57 public: 58 // Construction and destruction sparse_vector_element(vector_type & v,size_type i)59 sparse_vector_element (vector_type &v, size_type i): 60 container_reference<vector_type> (v), i_ (i) { 61 } 62 BOOST_UBLAS_INLINE sparse_vector_element(const sparse_vector_element & p)63 sparse_vector_element (const sparse_vector_element &p): 64 container_reference<vector_type> (p), i_ (p.i_) {} 65 BOOST_UBLAS_INLINE ~sparse_vector_element()66 ~sparse_vector_element () { 67 } 68 69 // Assignment 70 BOOST_UBLAS_INLINE operator =(const sparse_vector_element & p)71 sparse_vector_element &operator = (const sparse_vector_element &p) { 72 // Overide the implict copy assignment 73 p.get_d (); 74 set (p.d_); 75 return *this; 76 } 77 template<class D> 78 BOOST_UBLAS_INLINE operator =(const D & d)79 sparse_vector_element &operator = (const D &d) { 80 set (d); 81 return *this; 82 } 83 template<class D> 84 BOOST_UBLAS_INLINE operator +=(const D & d)85 sparse_vector_element &operator += (const D &d) { 86 get_d (); 87 d_ += d; 88 set (d_); 89 return *this; 90 } 91 template<class D> 92 BOOST_UBLAS_INLINE operator -=(const D & d)93 sparse_vector_element &operator -= (const D &d) { 94 get_d (); 95 d_ -= d; 96 set (d_); 97 return *this; 98 } 99 template<class D> 100 BOOST_UBLAS_INLINE operator *=(const D & d)101 sparse_vector_element &operator *= (const D &d) { 102 get_d (); 103 d_ *= d; 104 set (d_); 105 return *this; 106 } 107 template<class D> 108 BOOST_UBLAS_INLINE operator /=(const D & d)109 sparse_vector_element &operator /= (const D &d) { 110 get_d (); 111 d_ /= d; 112 set (d_); 113 return *this; 114 } 115 116 // Comparison 117 template<class D> 118 BOOST_UBLAS_INLINE operator ==(const D & d) const119 bool operator == (const D &d) const { 120 get_d (); 121 return d_ == d; 122 } 123 template<class D> 124 BOOST_UBLAS_INLINE operator !=(const D & d) const125 bool operator != (const D &d) const { 126 get_d (); 127 return d_ != d; 128 } 129 130 // Conversion - weak link in proxy as d_ is not a perfect alias for the element 131 BOOST_UBLAS_INLINE operator const_reference() const132 operator const_reference () const { 133 get_d (); 134 return d_; 135 } 136 137 // Conversion to reference - may be invalidated 138 BOOST_UBLAS_INLINE ref() const139 value_type& ref () const { 140 const pointer p = (*this) ().find_element (i_); 141 if (!p) 142 return (*this) ().insert_element (i_, value_type/*zero*/()); 143 else 144 return *p; 145 } 146 147 private: 148 size_type i_; 149 mutable value_type d_; 150 }; 151 152 /* 153 * Generalise explicit reference access 154 */ 155 namespace detail { 156 template <class R> 157 struct element_reference { 158 typedef R& reference; get_referenceboost::numeric::ublas::detail::element_reference159 static reference get_reference (reference r) 160 { 161 return r; 162 } 163 }; 164 template <class V> 165 struct element_reference<sparse_vector_element<V> > { 166 typedef typename V::value_type& reference; get_referenceboost::numeric::ublas::detail::element_reference167 static reference get_reference (const sparse_vector_element<V>& sve) 168 { 169 return sve.ref (); 170 } 171 }; 172 } 173 template <class VER> ref(VER & ver)174 typename detail::element_reference<VER>::reference ref (VER& ver) { 175 return detail::element_reference<VER>::get_reference (ver); 176 } 177 template <class VER> ref(const VER & ver)178 typename detail::element_reference<VER>::reference ref (const VER& ver) { 179 return detail::element_reference<VER>::get_reference (ver); 180 } 181 182 183 template<class V> 184 struct type_traits<sparse_vector_element<V> > { 185 typedef typename V::value_type element_type; 186 typedef type_traits<sparse_vector_element<V> > self_type; 187 typedef typename type_traits<element_type>::value_type value_type; 188 typedef typename type_traits<element_type>::const_reference const_reference; 189 typedef sparse_vector_element<V> reference; 190 typedef typename type_traits<element_type>::real_type real_type; 191 typedef typename type_traits<element_type>::precision_type precision_type; 192 193 static const unsigned plus_complexity = type_traits<element_type>::plus_complexity; 194 static const unsigned multiplies_complexity = type_traits<element_type>::multiplies_complexity; 195 196 static 197 BOOST_UBLAS_INLINE realboost::numeric::ublas::type_traits198 real_type real (const_reference t) { 199 return type_traits<element_type>::real (t); 200 } 201 static 202 BOOST_UBLAS_INLINE imagboost::numeric::ublas::type_traits203 real_type imag (const_reference t) { 204 return type_traits<element_type>::imag (t); 205 } 206 static 207 BOOST_UBLAS_INLINE conjboost::numeric::ublas::type_traits208 value_type conj (const_reference t) { 209 return type_traits<element_type>::conj (t); 210 } 211 212 static 213 BOOST_UBLAS_INLINE type_absboost::numeric::ublas::type_traits214 real_type type_abs (const_reference t) { 215 return type_traits<element_type>::type_abs (t); 216 } 217 static 218 BOOST_UBLAS_INLINE type_sqrtboost::numeric::ublas::type_traits219 value_type type_sqrt (const_reference t) { 220 return type_traits<element_type>::type_sqrt (t); 221 } 222 223 static 224 BOOST_UBLAS_INLINE norm_1boost::numeric::ublas::type_traits225 real_type norm_1 (const_reference t) { 226 return type_traits<element_type>::norm_1 (t); 227 } 228 static 229 BOOST_UBLAS_INLINE norm_2boost::numeric::ublas::type_traits230 real_type norm_2 (const_reference t) { 231 return type_traits<element_type>::norm_2 (t); 232 } 233 static 234 BOOST_UBLAS_INLINE norm_infboost::numeric::ublas::type_traits235 real_type norm_inf (const_reference t) { 236 return type_traits<element_type>::norm_inf (t); 237 } 238 239 static 240 BOOST_UBLAS_INLINE equalsboost::numeric::ublas::type_traits241 bool equals (const_reference t1, const_reference t2) { 242 return type_traits<element_type>::equals (t1, t2); 243 } 244 }; 245 246 template<class V1, class T2> 247 struct promote_traits<sparse_vector_element<V1>, T2> { 248 typedef typename promote_traits<typename sparse_vector_element<V1>::value_type, T2>::promote_type promote_type; 249 }; 250 template<class T1, class V2> 251 struct promote_traits<T1, sparse_vector_element<V2> > { 252 typedef typename promote_traits<T1, typename sparse_vector_element<V2>::value_type>::promote_type promote_type; 253 }; 254 template<class V1, class V2> 255 struct promote_traits<sparse_vector_element<V1>, sparse_vector_element<V2> > { 256 typedef typename promote_traits<typename sparse_vector_element<V1>::value_type, 257 typename sparse_vector_element<V2>::value_type>::promote_type promote_type; 258 }; 259 260 #endif 261 262 263 /** \brief Index map based sparse vector 264 * 265 * A sparse vector of values of type T of variable size. The sparse storage type A can be 266 * \c std::map<size_t, T> or \c map_array<size_t, T>. This means that only non-zero elements 267 * are effectively stored. 268 * 269 * For a \f$n\f$-dimensional sparse vector, and 0 <= i < n the non-zero elements \f$v_i\f$ 270 * are mapped to consecutive elements of the associative container, i.e. for elements 271 * \f$k = v_{i_1}\f$ and \f$k + 1 = v_{i_2}\f$ of the container, holds \f$i_1 < i_2\f$. 272 * 273 * Supported parameters for the adapted array are \c map_array<std::size_t, T> and 274 * \c map_std<std::size_t, T>. The latter is equivalent to \c std::map<std::size_t, T>. 275 * 276 * \tparam T the type of object stored in the vector (like double, float, complex, etc...) 277 * \tparam A the type of Storage array 278 */ 279 template<class T, class A> 280 class mapped_vector: 281 public vector_container<mapped_vector<T, A> > { 282 283 typedef T &true_reference; 284 typedef T *pointer; 285 typedef const T *const_pointer; 286 typedef mapped_vector<T, A> self_type; 287 public: 288 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS 289 using vector_container<self_type>::operator (); 290 #endif 291 typedef typename A::size_type size_type; 292 typedef typename A::difference_type difference_type; 293 typedef T value_type; 294 typedef A array_type; 295 typedef const value_type &const_reference; 296 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 297 typedef typename detail::map_traits<A,T>::reference reference; 298 #else 299 typedef sparse_vector_element<self_type> reference; 300 #endif 301 typedef const vector_reference<const self_type> const_closure_type; 302 typedef vector_reference<self_type> closure_type; 303 typedef self_type vector_temporary_type; 304 typedef sparse_tag storage_category; 305 306 // Construction and destruction 307 BOOST_UBLAS_INLINE mapped_vector()308 mapped_vector (): 309 vector_container<self_type> (), 310 size_ (0), data_ () {} 311 BOOST_UBLAS_INLINE mapped_vector(size_type size,size_type non_zeros=0)312 mapped_vector (size_type size, size_type non_zeros = 0): 313 vector_container<self_type> (), 314 size_ (size), data_ () { 315 detail::map_reserve (data(), restrict_capacity (non_zeros)); 316 } 317 BOOST_UBLAS_INLINE mapped_vector(const mapped_vector & v)318 mapped_vector (const mapped_vector &v): 319 vector_container<self_type> (), 320 size_ (v.size_), data_ (v.data_) {} 321 template<class AE> 322 BOOST_UBLAS_INLINE mapped_vector(const vector_expression<AE> & ae,size_type non_zeros=0)323 mapped_vector (const vector_expression<AE> &ae, size_type non_zeros = 0): 324 vector_container<self_type> (), 325 size_ (ae ().size ()), data_ () { 326 detail::map_reserve (data(), restrict_capacity (non_zeros)); 327 vector_assign<scalar_assign> (*this, ae); 328 } 329 330 // Accessors 331 BOOST_UBLAS_INLINE size() const332 size_type size () const { 333 return size_; 334 } 335 BOOST_UBLAS_INLINE nnz_capacity() const336 size_type nnz_capacity () const { 337 return detail::map_capacity (data ()); 338 } 339 BOOST_UBLAS_INLINE nnz() const340 size_type nnz () const { 341 return data (). size (); 342 } 343 344 // Storage accessors 345 BOOST_UBLAS_INLINE data() const346 const array_type &data () const { 347 return data_; 348 } 349 BOOST_UBLAS_INLINE data()350 array_type &data () { 351 return data_; 352 } 353 354 // Resizing 355 private: 356 BOOST_UBLAS_INLINE restrict_capacity(size_type non_zeros) const357 size_type restrict_capacity (size_type non_zeros) const { 358 non_zeros = (std::min) (non_zeros, size_); 359 return non_zeros; 360 } 361 public: 362 BOOST_UBLAS_INLINE resize(size_type size,bool preserve=true)363 void resize (size_type size, bool preserve = true) { 364 size_ = size; 365 if (preserve) { 366 data ().erase (data ().lower_bound(size_), data ().end()); 367 } 368 else { 369 data ().clear (); 370 } 371 } 372 373 // Reserving 374 BOOST_UBLAS_INLINE reserve(size_type non_zeros=0,bool preserve=true)375 void reserve (size_type non_zeros = 0, bool preserve = true) { 376 detail::map_reserve (data (), restrict_capacity (non_zeros)); 377 } 378 379 // Element support 380 BOOST_UBLAS_INLINE find_element(size_type i)381 pointer find_element (size_type i) { 382 return const_cast<pointer> (const_cast<const self_type&>(*this).find_element (i)); 383 } 384 BOOST_UBLAS_INLINE find_element(size_type i) const385 const_pointer find_element (size_type i) const { 386 const_subiterator_type it (data ().find (i)); 387 if (it == data ().end ()) 388 return 0; 389 BOOST_UBLAS_CHECK ((*it).first == i, internal_logic ()); // broken map 390 return &(*it).second; 391 } 392 393 // Element access 394 BOOST_UBLAS_INLINE operator ()(size_type i) const395 const_reference operator () (size_type i) const { 396 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 397 const_subiterator_type it (data ().find (i)); 398 if (it == data ().end ()) 399 return zero_; 400 BOOST_UBLAS_CHECK ((*it).first == i, internal_logic ()); // broken map 401 return (*it).second; 402 } 403 BOOST_UBLAS_INLINE ref(size_type i)404 true_reference ref (size_type i) { 405 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 406 std::pair<subiterator_type, bool> ii (data ().insert (typename array_type::value_type (i, value_type/*zero*/()))); 407 BOOST_UBLAS_CHECK ((ii.first)->first == i, internal_logic ()); // broken map 408 return (ii.first)->second; 409 } 410 BOOST_UBLAS_INLINE operator ()(size_type i)411 reference operator () (size_type i) { 412 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 413 return ref (i); 414 #else 415 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 416 return reference (*this, i); 417 #endif 418 } 419 420 BOOST_UBLAS_INLINE operator [](size_type i) const421 const_reference operator [] (size_type i) const { 422 return (*this) (i); 423 } 424 BOOST_UBLAS_INLINE operator [](size_type i)425 reference operator [] (size_type i) { 426 return (*this) (i); 427 } 428 429 // Element assignment 430 BOOST_UBLAS_INLINE insert_element(size_type i,const_reference t)431 true_reference insert_element (size_type i, const_reference t) { 432 std::pair<subiterator_type, bool> ii = data ().insert (typename array_type::value_type (i, t)); 433 BOOST_UBLAS_CHECK (ii.second, bad_index ()); // duplicate element 434 BOOST_UBLAS_CHECK ((ii.first)->first == i, internal_logic ()); // broken map 435 if (!ii.second) // existing element 436 (ii.first)->second = t; 437 return (ii.first)->second; 438 } 439 BOOST_UBLAS_INLINE erase_element(size_type i)440 void erase_element (size_type i) { 441 subiterator_type it = data ().find (i); 442 if (it == data ().end ()) 443 return; 444 data ().erase (it); 445 } 446 447 // Zeroing 448 BOOST_UBLAS_INLINE clear()449 void clear () { 450 data ().clear (); 451 } 452 453 // Assignment 454 BOOST_UBLAS_INLINE operator =(const mapped_vector & v)455 mapped_vector &operator = (const mapped_vector &v) { 456 if (this != &v) { 457 size_ = v.size_; 458 data () = v.data (); 459 } 460 return *this; 461 } 462 template<class C> // Container assignment without temporary 463 BOOST_UBLAS_INLINE operator =(const vector_container<C> & v)464 mapped_vector &operator = (const vector_container<C> &v) { 465 resize (v ().size (), false); 466 assign (v); 467 return *this; 468 } 469 BOOST_UBLAS_INLINE assign_temporary(mapped_vector & v)470 mapped_vector &assign_temporary (mapped_vector &v) { 471 swap (v); 472 return *this; 473 } 474 template<class AE> 475 BOOST_UBLAS_INLINE operator =(const vector_expression<AE> & ae)476 mapped_vector &operator = (const vector_expression<AE> &ae) { 477 self_type temporary (ae, detail::map_capacity (data())); 478 return assign_temporary (temporary); 479 } 480 template<class AE> 481 BOOST_UBLAS_INLINE assign(const vector_expression<AE> & ae)482 mapped_vector &assign (const vector_expression<AE> &ae) { 483 vector_assign<scalar_assign> (*this, ae); 484 return *this; 485 } 486 487 // Computed assignment 488 template<class AE> 489 BOOST_UBLAS_INLINE operator +=(const vector_expression<AE> & ae)490 mapped_vector &operator += (const vector_expression<AE> &ae) { 491 self_type temporary (*this + ae, detail::map_capacity (data())); 492 return assign_temporary (temporary); 493 } 494 template<class C> // Container assignment without temporary 495 BOOST_UBLAS_INLINE operator +=(const vector_container<C> & v)496 mapped_vector &operator += (const vector_container<C> &v) { 497 plus_assign (v); 498 return *this; 499 } 500 template<class AE> 501 BOOST_UBLAS_INLINE plus_assign(const vector_expression<AE> & ae)502 mapped_vector &plus_assign (const vector_expression<AE> &ae) { 503 vector_assign<scalar_plus_assign> (*this, ae); 504 return *this; 505 } 506 template<class AE> 507 BOOST_UBLAS_INLINE operator -=(const vector_expression<AE> & ae)508 mapped_vector &operator -= (const vector_expression<AE> &ae) { 509 self_type temporary (*this - ae, detail::map_capacity (data())); 510 return assign_temporary (temporary); 511 } 512 template<class C> // Container assignment without temporary 513 BOOST_UBLAS_INLINE operator -=(const vector_container<C> & v)514 mapped_vector &operator -= (const vector_container<C> &v) { 515 minus_assign (v); 516 return *this; 517 } 518 template<class AE> 519 BOOST_UBLAS_INLINE minus_assign(const vector_expression<AE> & ae)520 mapped_vector &minus_assign (const vector_expression<AE> &ae) { 521 vector_assign<scalar_minus_assign> (*this, ae); 522 return *this; 523 } 524 template<class AT> 525 BOOST_UBLAS_INLINE operator *=(const AT & at)526 mapped_vector &operator *= (const AT &at) { 527 vector_assign_scalar<scalar_multiplies_assign> (*this, at); 528 return *this; 529 } 530 template<class AT> 531 BOOST_UBLAS_INLINE operator /=(const AT & at)532 mapped_vector &operator /= (const AT &at) { 533 vector_assign_scalar<scalar_divides_assign> (*this, at); 534 return *this; 535 } 536 537 // Swapping 538 BOOST_UBLAS_INLINE swap(mapped_vector & v)539 void swap (mapped_vector &v) { 540 if (this != &v) { 541 std::swap (size_, v.size_); 542 data ().swap (v.data ()); 543 } 544 } 545 BOOST_UBLAS_INLINE swap(mapped_vector & v1,mapped_vector & v2)546 friend void swap (mapped_vector &v1, mapped_vector &v2) { 547 v1.swap (v2); 548 } 549 550 // Iterator types 551 private: 552 // Use storage iterator 553 typedef typename A::const_iterator const_subiterator_type; 554 typedef typename A::iterator subiterator_type; 555 556 BOOST_UBLAS_INLINE at_element(size_type i)557 true_reference at_element (size_type i) { 558 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 559 subiterator_type it (data ().find (i)); 560 BOOST_UBLAS_CHECK (it != data ().end(), bad_index ()); 561 BOOST_UBLAS_CHECK ((*it).first == i, internal_logic ()); // broken map 562 return it->second; 563 } 564 565 public: 566 class const_iterator; 567 class iterator; 568 569 // Element lookup 570 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i) const571 const_iterator find (size_type i) const { 572 return const_iterator (*this, data ().lower_bound (i)); 573 } 574 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i)575 iterator find (size_type i) { 576 return iterator (*this, data ().lower_bound (i)); 577 } 578 579 580 class const_iterator: 581 public container_const_reference<mapped_vector>, 582 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 583 const_iterator, value_type> { 584 public: 585 typedef typename mapped_vector::value_type value_type; 586 typedef typename mapped_vector::difference_type difference_type; 587 typedef typename mapped_vector::const_reference reference; 588 typedef const typename mapped_vector::pointer pointer; 589 590 // Construction and destruction 591 BOOST_UBLAS_INLINE const_iterator()592 const_iterator (): 593 container_const_reference<self_type> (), it_ () {} 594 BOOST_UBLAS_INLINE const_iterator(const self_type & v,const const_subiterator_type & it)595 const_iterator (const self_type &v, const const_subiterator_type &it): 596 container_const_reference<self_type> (v), it_ (it) {} 597 BOOST_UBLAS_INLINE const_iterator(const typename self_type::iterator & it)598 const_iterator (const typename self_type::iterator &it): // ISSUE self_type:: stops VC8 using std::iterator here 599 container_const_reference<self_type> (it ()), it_ (it.it_) {} 600 601 // Arithmetic 602 BOOST_UBLAS_INLINE operator ++()603 const_iterator &operator ++ () { 604 ++ it_; 605 return *this; 606 } 607 BOOST_UBLAS_INLINE operator --()608 const_iterator &operator -- () { 609 -- it_; 610 return *this; 611 } 612 613 // Dereference 614 BOOST_UBLAS_INLINE operator *() const615 const_reference operator * () const { 616 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 617 return (*it_).second; 618 } 619 620 // Index 621 BOOST_UBLAS_INLINE index() const622 size_type index () const { 623 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 624 BOOST_UBLAS_CHECK ((*it_).first < (*this) ().size (), bad_index ()); 625 return (*it_).first; 626 } 627 628 // Assignment 629 BOOST_UBLAS_INLINE operator =(const const_iterator & it)630 const_iterator &operator = (const const_iterator &it) { 631 container_const_reference<self_type>::assign (&it ()); 632 it_ = it.it_; 633 return *this; 634 } 635 636 // Comparison 637 BOOST_UBLAS_INLINE operator ==(const const_iterator & it) const638 bool operator == (const const_iterator &it) const { 639 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 640 return it_ == it.it_; 641 } 642 643 private: 644 const_subiterator_type it_; 645 }; 646 647 BOOST_UBLAS_INLINE begin() const648 const_iterator begin () const { 649 return const_iterator (*this, data ().begin ()); 650 } 651 BOOST_UBLAS_INLINE cbegin() const652 const_iterator cbegin () const { 653 return begin (); 654 } 655 BOOST_UBLAS_INLINE end() const656 const_iterator end () const { 657 return const_iterator (*this, data ().end ()); 658 } 659 BOOST_UBLAS_INLINE cend() const660 const_iterator cend () const { 661 return end (); 662 } 663 664 class iterator: 665 public container_reference<mapped_vector>, 666 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 667 iterator, value_type> { 668 public: 669 typedef typename mapped_vector::value_type value_type; 670 typedef typename mapped_vector::difference_type difference_type; 671 typedef typename mapped_vector::true_reference reference; 672 typedef typename mapped_vector::pointer pointer; 673 674 // Construction and destruction 675 BOOST_UBLAS_INLINE iterator()676 iterator (): 677 container_reference<self_type> (), it_ () {} 678 BOOST_UBLAS_INLINE iterator(self_type & v,const subiterator_type & it)679 iterator (self_type &v, const subiterator_type &it): 680 container_reference<self_type> (v), it_ (it) {} 681 682 // Arithmetic 683 BOOST_UBLAS_INLINE operator ++()684 iterator &operator ++ () { 685 ++ it_; 686 return *this; 687 } 688 BOOST_UBLAS_INLINE operator --()689 iterator &operator -- () { 690 -- it_; 691 return *this; 692 } 693 694 // Dereference 695 BOOST_UBLAS_INLINE operator *() const696 reference operator * () const { 697 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 698 return (*it_).second; 699 } 700 701 // Index 702 BOOST_UBLAS_INLINE index() const703 size_type index () const { 704 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 705 BOOST_UBLAS_CHECK ((*it_).first < (*this) ().size (), bad_index ()); 706 return (*it_).first; 707 } 708 709 // Assignment 710 BOOST_UBLAS_INLINE operator =(const iterator & it)711 iterator &operator = (const iterator &it) { 712 container_reference<self_type>::assign (&it ()); 713 it_ = it.it_; 714 return *this; 715 } 716 717 // Comparison 718 BOOST_UBLAS_INLINE operator ==(const iterator & it) const719 bool operator == (const iterator &it) const { 720 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 721 return it_ == it.it_; 722 } 723 724 private: 725 subiterator_type it_; 726 727 friend class const_iterator; 728 }; 729 730 BOOST_UBLAS_INLINE begin()731 iterator begin () { 732 return iterator (*this, data ().begin ()); 733 } 734 BOOST_UBLAS_INLINE end()735 iterator end () { 736 return iterator (*this, data ().end ()); 737 } 738 739 // Reverse iterator 740 typedef reverse_iterator_base<const_iterator> const_reverse_iterator; 741 typedef reverse_iterator_base<iterator> reverse_iterator; 742 743 BOOST_UBLAS_INLINE rbegin() const744 const_reverse_iterator rbegin () const { 745 return const_reverse_iterator (end ()); 746 } 747 BOOST_UBLAS_INLINE crbegin() const748 const_reverse_iterator crbegin () const { 749 return rbegin (); 750 } 751 BOOST_UBLAS_INLINE rend() const752 const_reverse_iterator rend () const { 753 return const_reverse_iterator (begin ()); 754 } 755 BOOST_UBLAS_INLINE crend() const756 const_reverse_iterator crend () const { 757 return rend (); 758 } 759 BOOST_UBLAS_INLINE rbegin()760 reverse_iterator rbegin () { 761 return reverse_iterator (end ()); 762 } 763 BOOST_UBLAS_INLINE rend()764 reverse_iterator rend () { 765 return reverse_iterator (begin ()); 766 } 767 768 // Serialization 769 template<class Archive> serialize(Archive & ar,const unsigned int)770 void serialize(Archive & ar, const unsigned int /* file_version */){ 771 serialization::collection_size_type s (size_); 772 ar & serialization::make_nvp("size",s); 773 if (Archive::is_loading::value) { 774 size_ = s; 775 } 776 ar & serialization::make_nvp("data", data_); 777 } 778 779 private: 780 size_type size_; 781 array_type data_; 782 static const value_type zero_; 783 }; 784 785 template<class T, class A> 786 const typename mapped_vector<T, A>::value_type mapped_vector<T, A>::zero_ = value_type/*zero*/(); 787 788 789 // Thanks to Kresimir Fresl for extending this to cover different index bases. 790 791 /** \brief Compressed array based sparse vector 792 * 793 * a sparse vector of values of type T of variable size. The non zero values are stored as 794 * two seperate arrays: an index array and a value array. The index array is always sorted 795 * and there is at most one entry for each index. Inserting an element can be time consuming. 796 * If the vector contains a few zero entries, then it is better to have a normal vector. 797 * If the vector has a very high dimension with a few non-zero values, then this vector is 798 * very memory efficient (at the cost of a few more computations). 799 * 800 * For a \f$n\f$-dimensional compressed vector and \f$0 \leq i < n\f$ the non-zero elements 801 * \f$v_i\f$ are mapped to consecutive elements of the index and value container, i.e. for 802 * elements \f$k = v_{i_1}\f$ and \f$k + 1 = v_{i_2}\f$ of these containers holds \f$i_1 < i_2\f$. 803 * 804 * Supported parameters for the adapted array (indices and values) are \c unbounded_array<> , 805 * \c bounded_array<> and \c std::vector<>. 806 * 807 * \tparam T the type of object stored in the vector (like double, float, complex, etc...) 808 * \tparam IB the index base of the compressed vector. Default is 0. Other supported value is 1 809 * \tparam IA the type of adapted array for indices. Default is \c unbounded_array<std::size_t> 810 * \tparam TA the type of adapted array for values. Default is unbounded_array<T> 811 */ 812 template<class T, std::size_t IB, class IA, class TA> 813 class compressed_vector: 814 public vector_container<compressed_vector<T, IB, IA, TA> > { 815 816 typedef T &true_reference; 817 typedef T *pointer; 818 typedef const T *const_pointer; 819 typedef compressed_vector<T, IB, IA, TA> self_type; 820 public: 821 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS 822 using vector_container<self_type>::operator (); 823 #endif 824 // ISSUE require type consistency check 825 // is_convertable (IA::size_type, TA::size_type) 826 typedef typename IA::value_type size_type; 827 typedef typename IA::difference_type difference_type; 828 typedef T value_type; 829 typedef const T &const_reference; 830 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 831 typedef T &reference; 832 #else 833 typedef sparse_vector_element<self_type> reference; 834 #endif 835 typedef IA index_array_type; 836 typedef TA value_array_type; 837 typedef const vector_reference<const self_type> const_closure_type; 838 typedef vector_reference<self_type> closure_type; 839 typedef self_type vector_temporary_type; 840 typedef sparse_tag storage_category; 841 842 // Construction and destruction 843 BOOST_UBLAS_INLINE compressed_vector()844 compressed_vector (): 845 vector_container<self_type> (), 846 size_ (0), capacity_ (restrict_capacity (0)), filled_ (0), 847 index_data_ (capacity_), value_data_ (capacity_) { 848 storage_invariants (); 849 } 850 explicit BOOST_UBLAS_INLINE compressed_vector(size_type size,size_type non_zeros=0)851 compressed_vector (size_type size, size_type non_zeros = 0): 852 vector_container<self_type> (), 853 size_ (size), capacity_ (restrict_capacity (non_zeros)), filled_ (0), 854 index_data_ (capacity_), value_data_ (capacity_) { 855 storage_invariants (); 856 } 857 BOOST_UBLAS_INLINE compressed_vector(const compressed_vector & v)858 compressed_vector (const compressed_vector &v): 859 vector_container<self_type> (), 860 size_ (v.size_), capacity_ (v.capacity_), filled_ (v.filled_), 861 index_data_ (v.index_data_), value_data_ (v.value_data_) { 862 storage_invariants (); 863 } 864 template<class AE> 865 BOOST_UBLAS_INLINE compressed_vector(const vector_expression<AE> & ae,size_type non_zeros=0)866 compressed_vector (const vector_expression<AE> &ae, size_type non_zeros = 0): 867 vector_container<self_type> (), 868 size_ (ae ().size ()), capacity_ (restrict_capacity (non_zeros)), filled_ (0), 869 index_data_ (capacity_), value_data_ (capacity_) { 870 storage_invariants (); 871 vector_assign<scalar_assign> (*this, ae); 872 } 873 874 // Accessors 875 BOOST_UBLAS_INLINE size() const876 size_type size () const { 877 return size_; 878 } 879 BOOST_UBLAS_INLINE nnz_capacity() const880 size_type nnz_capacity () const { 881 return capacity_; 882 } 883 BOOST_UBLAS_INLINE nnz() const884 size_type nnz () const { 885 return filled_; 886 } 887 888 // Storage accessors 889 BOOST_UBLAS_INLINE index_base()890 static size_type index_base () { 891 return IB; 892 } 893 BOOST_UBLAS_INLINE filled() const894 typename index_array_type::size_type filled () const { 895 return filled_; 896 } 897 BOOST_UBLAS_INLINE index_data() const898 const index_array_type &index_data () const { 899 return index_data_; 900 } 901 BOOST_UBLAS_INLINE value_data() const902 const value_array_type &value_data () const { 903 return value_data_; 904 } 905 BOOST_UBLAS_INLINE set_filled(const typename index_array_type::size_type & filled)906 void set_filled (const typename index_array_type::size_type & filled) { 907 filled_ = filled; 908 storage_invariants (); 909 } 910 BOOST_UBLAS_INLINE index_data()911 index_array_type &index_data () { 912 return index_data_; 913 } 914 BOOST_UBLAS_INLINE value_data()915 value_array_type &value_data () { 916 return value_data_; 917 } 918 919 // Resizing 920 private: 921 BOOST_UBLAS_INLINE restrict_capacity(size_type non_zeros) const922 size_type restrict_capacity (size_type non_zeros) const { 923 non_zeros = (std::max) (non_zeros, size_type (1)); 924 non_zeros = (std::min) (non_zeros, size_); 925 return non_zeros; 926 } 927 public: 928 BOOST_UBLAS_INLINE resize(size_type size,bool preserve=true)929 void resize (size_type size, bool preserve = true) { 930 size_ = size; 931 capacity_ = restrict_capacity (capacity_); 932 if (preserve) { 933 index_data_. resize (capacity_, size_type ()); 934 value_data_. resize (capacity_, value_type ()); 935 filled_ = (std::min) (capacity_, filled_); 936 while ((filled_ > 0) && (zero_based(index_data_[filled_ - 1]) >= size)) { 937 --filled_; 938 } 939 } 940 else { 941 index_data_. resize (capacity_); 942 value_data_. resize (capacity_); 943 filled_ = 0; 944 } 945 storage_invariants (); 946 } 947 948 // Reserving 949 BOOST_UBLAS_INLINE reserve(size_type non_zeros,bool preserve=true)950 void reserve (size_type non_zeros, bool preserve = true) { 951 capacity_ = restrict_capacity (non_zeros); 952 if (preserve) { 953 index_data_. resize (capacity_, size_type ()); 954 value_data_. resize (capacity_, value_type ()); 955 filled_ = (std::min) (capacity_, filled_); 956 } 957 else { 958 index_data_. resize (capacity_); 959 value_data_. resize (capacity_); 960 filled_ = 0; 961 } 962 storage_invariants (); 963 } 964 965 // Element support 966 BOOST_UBLAS_INLINE find_element(size_type i)967 pointer find_element (size_type i) { 968 return const_cast<pointer> (const_cast<const self_type&>(*this).find_element (i)); 969 } 970 BOOST_UBLAS_INLINE find_element(size_type i) const971 const_pointer find_element (size_type i) const { 972 const_subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 973 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 974 return 0; 975 return &value_data_ [it - index_data_.begin ()]; 976 } 977 978 // Element access 979 BOOST_UBLAS_INLINE operator ()(size_type i) const980 const_reference operator () (size_type i) const { 981 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 982 const_subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 983 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 984 return zero_; 985 return value_data_ [it - index_data_.begin ()]; 986 } 987 BOOST_UBLAS_INLINE ref(size_type i)988 true_reference ref (size_type i) { 989 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 990 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 991 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 992 return insert_element (i, value_type/*zero*/()); 993 else 994 return value_data_ [it - index_data_.begin ()]; 995 } 996 BOOST_UBLAS_INLINE operator ()(size_type i)997 reference operator () (size_type i) { 998 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 999 return ref (i) ; 1000 #else 1001 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1002 return reference (*this, i); 1003 #endif 1004 } 1005 1006 BOOST_UBLAS_INLINE operator [](size_type i) const1007 const_reference operator [] (size_type i) const { 1008 return (*this) (i); 1009 } 1010 BOOST_UBLAS_INLINE operator [](size_type i)1011 reference operator [] (size_type i) { 1012 return (*this) (i); 1013 } 1014 1015 // Element assignment 1016 BOOST_UBLAS_INLINE insert_element(size_type i,const_reference t)1017 true_reference insert_element (size_type i, const_reference t) { 1018 BOOST_UBLAS_CHECK (!find_element (i), bad_index ()); // duplicate element 1019 if (filled_ >= capacity_) 1020 reserve (2 * capacity_, true); 1021 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1022 // ISSUE max_capacity limit due to difference_type 1023 typename std::iterator_traits<subiterator_type>::difference_type n = it - index_data_.begin (); 1024 BOOST_UBLAS_CHECK (filled_ == 0 || filled_ == typename index_array_type::size_type (n) || *it != k_based (i), internal_logic ()); // duplicate found by lower_bound 1025 ++ filled_; 1026 it = index_data_.begin () + n; 1027 std::copy_backward (it, index_data_.begin () + filled_ - 1, index_data_.begin () + filled_); 1028 *it = k_based (i); 1029 typename value_array_type::iterator itt (value_data_.begin () + n); 1030 std::copy_backward (itt, value_data_.begin () + filled_ - 1, value_data_.begin () + filled_); 1031 *itt = t; 1032 storage_invariants (); 1033 return *itt; 1034 } 1035 BOOST_UBLAS_INLINE erase_element(size_type i)1036 void erase_element (size_type i) { 1037 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1038 typename std::iterator_traits<subiterator_type>::difference_type n = it - index_data_.begin (); 1039 if (filled_ > typename index_array_type::size_type (n) && *it == k_based (i)) { 1040 std::copy (it + 1, index_data_.begin () + filled_, it); 1041 typename value_array_type::iterator itt (value_data_.begin () + n); 1042 std::copy (itt + 1, value_data_.begin () + filled_, itt); 1043 -- filled_; 1044 } 1045 storage_invariants (); 1046 } 1047 1048 // Zeroing 1049 BOOST_UBLAS_INLINE clear()1050 void clear () { 1051 filled_ = 0; 1052 storage_invariants (); 1053 } 1054 1055 // Assignment 1056 BOOST_UBLAS_INLINE operator =(const compressed_vector & v)1057 compressed_vector &operator = (const compressed_vector &v) { 1058 if (this != &v) { 1059 size_ = v.size_; 1060 capacity_ = v.capacity_; 1061 filled_ = v.filled_; 1062 index_data_ = v.index_data_; 1063 value_data_ = v.value_data_; 1064 } 1065 storage_invariants (); 1066 return *this; 1067 } 1068 template<class C> // Container assignment without temporary 1069 BOOST_UBLAS_INLINE operator =(const vector_container<C> & v)1070 compressed_vector &operator = (const vector_container<C> &v) { 1071 resize (v ().size (), false); 1072 assign (v); 1073 return *this; 1074 } 1075 BOOST_UBLAS_INLINE assign_temporary(compressed_vector & v)1076 compressed_vector &assign_temporary (compressed_vector &v) { 1077 swap (v); 1078 return *this; 1079 } 1080 template<class AE> 1081 BOOST_UBLAS_INLINE operator =(const vector_expression<AE> & ae)1082 compressed_vector &operator = (const vector_expression<AE> &ae) { 1083 self_type temporary (ae, capacity_); 1084 return assign_temporary (temporary); 1085 } 1086 template<class AE> 1087 BOOST_UBLAS_INLINE assign(const vector_expression<AE> & ae)1088 compressed_vector &assign (const vector_expression<AE> &ae) { 1089 vector_assign<scalar_assign> (*this, ae); 1090 return *this; 1091 } 1092 1093 // Computed assignment 1094 template<class AE> 1095 BOOST_UBLAS_INLINE operator +=(const vector_expression<AE> & ae)1096 compressed_vector &operator += (const vector_expression<AE> &ae) { 1097 self_type temporary (*this + ae, capacity_); 1098 return assign_temporary (temporary); 1099 } 1100 template<class C> // Container assignment without temporary 1101 BOOST_UBLAS_INLINE operator +=(const vector_container<C> & v)1102 compressed_vector &operator += (const vector_container<C> &v) { 1103 plus_assign (v); 1104 return *this; 1105 } 1106 template<class AE> 1107 BOOST_UBLAS_INLINE plus_assign(const vector_expression<AE> & ae)1108 compressed_vector &plus_assign (const vector_expression<AE> &ae) { 1109 vector_assign<scalar_plus_assign> (*this, ae); 1110 return *this; 1111 } 1112 template<class AE> 1113 BOOST_UBLAS_INLINE operator -=(const vector_expression<AE> & ae)1114 compressed_vector &operator -= (const vector_expression<AE> &ae) { 1115 self_type temporary (*this - ae, capacity_); 1116 return assign_temporary (temporary); 1117 } 1118 template<class C> // Container assignment without temporary 1119 BOOST_UBLAS_INLINE operator -=(const vector_container<C> & v)1120 compressed_vector &operator -= (const vector_container<C> &v) { 1121 minus_assign (v); 1122 return *this; 1123 } 1124 template<class AE> 1125 BOOST_UBLAS_INLINE minus_assign(const vector_expression<AE> & ae)1126 compressed_vector &minus_assign (const vector_expression<AE> &ae) { 1127 vector_assign<scalar_minus_assign> (*this, ae); 1128 return *this; 1129 } 1130 template<class AT> 1131 BOOST_UBLAS_INLINE operator *=(const AT & at)1132 compressed_vector &operator *= (const AT &at) { 1133 vector_assign_scalar<scalar_multiplies_assign> (*this, at); 1134 return *this; 1135 } 1136 template<class AT> 1137 BOOST_UBLAS_INLINE operator /=(const AT & at)1138 compressed_vector &operator /= (const AT &at) { 1139 vector_assign_scalar<scalar_divides_assign> (*this, at); 1140 return *this; 1141 } 1142 1143 // Swapping 1144 BOOST_UBLAS_INLINE swap(compressed_vector & v)1145 void swap (compressed_vector &v) { 1146 if (this != &v) { 1147 std::swap (size_, v.size_); 1148 std::swap (capacity_, v.capacity_); 1149 std::swap (filled_, v.filled_); 1150 index_data_.swap (v.index_data_); 1151 value_data_.swap (v.value_data_); 1152 } 1153 storage_invariants (); 1154 } 1155 BOOST_UBLAS_INLINE swap(compressed_vector & v1,compressed_vector & v2)1156 friend void swap (compressed_vector &v1, compressed_vector &v2) { 1157 v1.swap (v2); 1158 } 1159 1160 // Back element insertion and erasure 1161 BOOST_UBLAS_INLINE push_back(size_type i,const_reference t)1162 void push_back (size_type i, const_reference t) { 1163 BOOST_UBLAS_CHECK (filled_ == 0 || index_data_ [filled_ - 1] < k_based (i), external_logic ()); 1164 if (filled_ >= capacity_) 1165 reserve (2 * capacity_, true); 1166 BOOST_UBLAS_CHECK (filled_ < capacity_, internal_logic ()); 1167 index_data_ [filled_] = k_based (i); 1168 value_data_ [filled_] = t; 1169 ++ filled_; 1170 storage_invariants (); 1171 } 1172 BOOST_UBLAS_INLINE pop_back()1173 void pop_back () { 1174 BOOST_UBLAS_CHECK (filled_ > 0, external_logic ()); 1175 -- filled_; 1176 storage_invariants (); 1177 } 1178 1179 // Iterator types 1180 private: 1181 // Use index array iterator 1182 typedef typename IA::const_iterator const_subiterator_type; 1183 typedef typename IA::iterator subiterator_type; 1184 1185 BOOST_UBLAS_INLINE at_element(size_type i)1186 true_reference at_element (size_type i) { 1187 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1188 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1189 BOOST_UBLAS_CHECK (it != index_data_.begin () + filled_ && *it == k_based (i), bad_index ()); 1190 return value_data_ [it - index_data_.begin ()]; 1191 } 1192 1193 public: 1194 class const_iterator; 1195 class iterator; 1196 1197 // Element lookup 1198 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i) const1199 const_iterator find (size_type i) const { 1200 return const_iterator (*this, detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1201 } 1202 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i)1203 iterator find (size_type i) { 1204 return iterator (*this, detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1205 } 1206 1207 1208 class const_iterator: 1209 public container_const_reference<compressed_vector>, 1210 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 1211 const_iterator, value_type> { 1212 public: 1213 typedef typename compressed_vector::value_type value_type; 1214 typedef typename compressed_vector::difference_type difference_type; 1215 typedef typename compressed_vector::const_reference reference; 1216 typedef const typename compressed_vector::pointer pointer; 1217 1218 // Construction and destruction 1219 BOOST_UBLAS_INLINE const_iterator()1220 const_iterator (): 1221 container_const_reference<self_type> (), it_ () {} 1222 BOOST_UBLAS_INLINE const_iterator(const self_type & v,const const_subiterator_type & it)1223 const_iterator (const self_type &v, const const_subiterator_type &it): 1224 container_const_reference<self_type> (v), it_ (it) {} 1225 BOOST_UBLAS_INLINE const_iterator(const typename self_type::iterator & it)1226 const_iterator (const typename self_type::iterator &it): // ISSUE self_type:: stops VC8 using std::iterator here 1227 container_const_reference<self_type> (it ()), it_ (it.it_) {} 1228 1229 // Arithmetic 1230 BOOST_UBLAS_INLINE operator ++()1231 const_iterator &operator ++ () { 1232 ++ it_; 1233 return *this; 1234 } 1235 BOOST_UBLAS_INLINE operator --()1236 const_iterator &operator -- () { 1237 -- it_; 1238 return *this; 1239 } 1240 1241 // Dereference 1242 BOOST_UBLAS_INLINE operator *() const1243 const_reference operator * () const { 1244 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 1245 return (*this) ().value_data_ [it_ - (*this) ().index_data_.begin ()]; 1246 } 1247 1248 // Index 1249 BOOST_UBLAS_INLINE index() const1250 size_type index () const { 1251 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 1252 BOOST_UBLAS_CHECK ((*this) ().zero_based (*it_) < (*this) ().size (), bad_index ()); 1253 return (*this) ().zero_based (*it_); 1254 } 1255 1256 // Assignment 1257 BOOST_UBLAS_INLINE operator =(const const_iterator & it)1258 const_iterator &operator = (const const_iterator &it) { 1259 container_const_reference<self_type>::assign (&it ()); 1260 it_ = it.it_; 1261 return *this; 1262 } 1263 1264 // Comparison 1265 BOOST_UBLAS_INLINE operator ==(const const_iterator & it) const1266 bool operator == (const const_iterator &it) const { 1267 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 1268 return it_ == it.it_; 1269 } 1270 1271 private: 1272 const_subiterator_type it_; 1273 }; 1274 1275 BOOST_UBLAS_INLINE begin() const1276 const_iterator begin () const { 1277 return find (0); 1278 } 1279 BOOST_UBLAS_INLINE cbegin() const1280 const_iterator cbegin () const { 1281 return begin (); 1282 } 1283 BOOST_UBLAS_INLINE end() const1284 const_iterator end () const { 1285 return find (size_); 1286 } 1287 BOOST_UBLAS_INLINE cend() const1288 const_iterator cend () const { 1289 return end (); 1290 } 1291 1292 class iterator: 1293 public container_reference<compressed_vector>, 1294 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 1295 iterator, value_type> { 1296 public: 1297 typedef typename compressed_vector::value_type value_type; 1298 typedef typename compressed_vector::difference_type difference_type; 1299 typedef typename compressed_vector::true_reference reference; 1300 typedef typename compressed_vector::pointer pointer; 1301 1302 // Construction and destruction 1303 BOOST_UBLAS_INLINE iterator()1304 iterator (): 1305 container_reference<self_type> (), it_ () {} 1306 BOOST_UBLAS_INLINE iterator(self_type & v,const subiterator_type & it)1307 iterator (self_type &v, const subiterator_type &it): 1308 container_reference<self_type> (v), it_ (it) {} 1309 1310 // Arithmetic 1311 BOOST_UBLAS_INLINE operator ++()1312 iterator &operator ++ () { 1313 ++ it_; 1314 return *this; 1315 } 1316 BOOST_UBLAS_INLINE operator --()1317 iterator &operator -- () { 1318 -- it_; 1319 return *this; 1320 } 1321 1322 // Dereference 1323 BOOST_UBLAS_INLINE operator *() const1324 reference operator * () const { 1325 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 1326 return (*this) ().value_data_ [it_ - (*this) ().index_data_.begin ()]; 1327 } 1328 1329 // Index 1330 BOOST_UBLAS_INLINE index() const1331 size_type index () const { 1332 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 1333 BOOST_UBLAS_CHECK ((*this) ().zero_based (*it_) < (*this) ().size (), bad_index ()); 1334 return (*this) ().zero_based (*it_); 1335 } 1336 1337 // Assignment 1338 BOOST_UBLAS_INLINE operator =(const iterator & it)1339 iterator &operator = (const iterator &it) { 1340 container_reference<self_type>::assign (&it ()); 1341 it_ = it.it_; 1342 return *this; 1343 } 1344 1345 // Comparison 1346 BOOST_UBLAS_INLINE operator ==(const iterator & it) const1347 bool operator == (const iterator &it) const { 1348 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 1349 return it_ == it.it_; 1350 } 1351 1352 private: 1353 subiterator_type it_; 1354 1355 friend class const_iterator; 1356 }; 1357 1358 BOOST_UBLAS_INLINE begin()1359 iterator begin () { 1360 return find (0); 1361 } 1362 BOOST_UBLAS_INLINE end()1363 iterator end () { 1364 return find (size_); 1365 } 1366 1367 // Reverse iterator 1368 typedef reverse_iterator_base<const_iterator> const_reverse_iterator; 1369 typedef reverse_iterator_base<iterator> reverse_iterator; 1370 1371 BOOST_UBLAS_INLINE rbegin() const1372 const_reverse_iterator rbegin () const { 1373 return const_reverse_iterator (end ()); 1374 } 1375 BOOST_UBLAS_INLINE crbegin() const1376 const_reverse_iterator crbegin () const { 1377 return rbegin (); 1378 } 1379 BOOST_UBLAS_INLINE rend() const1380 const_reverse_iterator rend () const { 1381 return const_reverse_iterator (begin ()); 1382 } 1383 BOOST_UBLAS_INLINE crend() const1384 const_reverse_iterator crend () const { 1385 return rend (); 1386 } 1387 BOOST_UBLAS_INLINE rbegin()1388 reverse_iterator rbegin () { 1389 return reverse_iterator (end ()); 1390 } 1391 BOOST_UBLAS_INLINE rend()1392 reverse_iterator rend () { 1393 return reverse_iterator (begin ()); 1394 } 1395 1396 // Serialization 1397 template<class Archive> serialize(Archive & ar,const unsigned int)1398 void serialize(Archive & ar, const unsigned int /* file_version */){ 1399 serialization::collection_size_type s (size_); 1400 ar & serialization::make_nvp("size",s); 1401 if (Archive::is_loading::value) { 1402 size_ = s; 1403 } 1404 // ISSUE: filled may be much less than capacity 1405 // ISSUE: index_data_ and value_data_ are undefined between filled and capacity (trouble with 'nan'-values) 1406 ar & serialization::make_nvp("capacity", capacity_); 1407 ar & serialization::make_nvp("filled", filled_); 1408 ar & serialization::make_nvp("index_data", index_data_); 1409 ar & serialization::make_nvp("value_data", value_data_); 1410 storage_invariants(); 1411 } 1412 1413 private: storage_invariants() const1414 void storage_invariants () const 1415 { 1416 BOOST_UBLAS_CHECK (capacity_ == index_data_.size (), internal_logic ()); 1417 BOOST_UBLAS_CHECK (capacity_ == value_data_.size (), internal_logic ()); 1418 BOOST_UBLAS_CHECK (filled_ <= capacity_, internal_logic ()); 1419 BOOST_UBLAS_CHECK ((0 == filled_) || (zero_based(index_data_[filled_ - 1]) < size_), internal_logic ()); 1420 } 1421 1422 size_type size_; 1423 typename index_array_type::size_type capacity_; 1424 typename index_array_type::size_type filled_; 1425 index_array_type index_data_; 1426 value_array_type value_data_; 1427 static const value_type zero_; 1428 1429 BOOST_UBLAS_INLINE zero_based(size_type k_based_index)1430 static size_type zero_based (size_type k_based_index) { 1431 return k_based_index - IB; 1432 } 1433 BOOST_UBLAS_INLINE k_based(size_type zero_based_index)1434 static size_type k_based (size_type zero_based_index) { 1435 return zero_based_index + IB; 1436 } 1437 1438 friend class iterator; 1439 friend class const_iterator; 1440 }; 1441 1442 template<class T, std::size_t IB, class IA, class TA> 1443 const typename compressed_vector<T, IB, IA, TA>::value_type compressed_vector<T, IB, IA, TA>::zero_ = value_type/*zero*/(); 1444 1445 // Thanks to Kresimir Fresl for extending this to cover different index bases. 1446 1447 /** \brief Coordimate array based sparse vector 1448 * 1449 * a sparse vector of values of type \c T of variable size. The non zero values are stored 1450 * as two seperate arrays: an index array and a value array. The arrays may be out of order 1451 * with multiple entries for each vector element. If there are multiple values for the same 1452 * index the sum of these values is the real value. It is way more efficient for inserting values 1453 * than a \c compressed_vector but less memory efficient. Also linearly parsing a vector can 1454 * be longer in specific cases than a \c compressed_vector. 1455 * 1456 * For a n-dimensional sorted coordinate vector and \f$ 0 \leq i < n\f$ the non-zero elements 1457 * \f$v_i\f$ are mapped to consecutive elements of the index and value container, i.e. for 1458 * elements \f$k = v_{i_1}\f$ and \f$k + 1 = v_{i_2}\f$ of these containers holds \f$i_1 < i_2\f$. 1459 * 1460 * Supported parameters for the adapted array (indices and values) are \c unbounded_array<> , 1461 * \c bounded_array<> and \c std::vector<>. 1462 * 1463 * \tparam T the type of object stored in the vector (like double, float, complex, etc...) 1464 * \tparam IB the index base of the compressed vector. Default is 0. Other supported value is 1 1465 * \tparam IA the type of adapted array for indices. Default is \c unbounded_array<std::size_t> 1466 * \tparam TA the type of adapted array for values. Default is unbounded_array<T> 1467 */ 1468 template<class T, std::size_t IB, class IA, class TA> 1469 class coordinate_vector: 1470 public vector_container<coordinate_vector<T, IB, IA, TA> > { 1471 1472 typedef T &true_reference; 1473 typedef T *pointer; 1474 typedef const T *const_pointer; 1475 typedef coordinate_vector<T, IB, IA, TA> self_type; 1476 public: 1477 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS 1478 using vector_container<self_type>::operator (); 1479 #endif 1480 // ISSUE require type consistency check 1481 // is_convertable (IA::size_type, TA::size_type) 1482 typedef typename IA::value_type size_type; 1483 typedef typename IA::difference_type difference_type; 1484 typedef T value_type; 1485 typedef const T &const_reference; 1486 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 1487 typedef T &reference; 1488 #else 1489 typedef sparse_vector_element<self_type> reference; 1490 #endif 1491 typedef IA index_array_type; 1492 typedef TA value_array_type; 1493 typedef const vector_reference<const self_type> const_closure_type; 1494 typedef vector_reference<self_type> closure_type; 1495 typedef self_type vector_temporary_type; 1496 typedef sparse_tag storage_category; 1497 1498 // Construction and destruction 1499 BOOST_UBLAS_INLINE coordinate_vector()1500 coordinate_vector (): 1501 vector_container<self_type> (), 1502 size_ (0), capacity_ (restrict_capacity (0)), 1503 filled_ (0), sorted_filled_ (filled_), sorted_ (true), 1504 index_data_ (capacity_), value_data_ (capacity_) { 1505 storage_invariants (); 1506 } 1507 explicit BOOST_UBLAS_INLINE coordinate_vector(size_type size,size_type non_zeros=0)1508 coordinate_vector (size_type size, size_type non_zeros = 0): 1509 vector_container<self_type> (), 1510 size_ (size), capacity_ (restrict_capacity (non_zeros)), 1511 filled_ (0), sorted_filled_ (filled_), sorted_ (true), 1512 index_data_ (capacity_), value_data_ (capacity_) { 1513 storage_invariants (); 1514 } 1515 BOOST_UBLAS_INLINE coordinate_vector(const coordinate_vector & v)1516 coordinate_vector (const coordinate_vector &v): 1517 vector_container<self_type> (), 1518 size_ (v.size_), capacity_ (v.capacity_), 1519 filled_ (v.filled_), sorted_filled_ (v.sorted_filled_), sorted_ (v.sorted_), 1520 index_data_ (v.index_data_), value_data_ (v.value_data_) { 1521 storage_invariants (); 1522 } 1523 template<class AE> 1524 BOOST_UBLAS_INLINE coordinate_vector(const vector_expression<AE> & ae,size_type non_zeros=0)1525 coordinate_vector (const vector_expression<AE> &ae, size_type non_zeros = 0): 1526 vector_container<self_type> (), 1527 size_ (ae ().size ()), capacity_ (restrict_capacity (non_zeros)), 1528 filled_ (0), sorted_filled_ (filled_), sorted_ (true), 1529 index_data_ (capacity_), value_data_ (capacity_) { 1530 storage_invariants (); 1531 vector_assign<scalar_assign> (*this, ae); 1532 } 1533 1534 // Accessors 1535 BOOST_UBLAS_INLINE size() const1536 size_type size () const { 1537 return size_; 1538 } 1539 BOOST_UBLAS_INLINE nnz_capacity() const1540 size_type nnz_capacity () const { 1541 return capacity_; 1542 } 1543 BOOST_UBLAS_INLINE nnz() const1544 size_type nnz () const { 1545 return filled_; 1546 } 1547 1548 // Storage accessors 1549 BOOST_UBLAS_INLINE index_base()1550 static size_type index_base () { 1551 return IB; 1552 } 1553 BOOST_UBLAS_INLINE filled() const1554 typename index_array_type::size_type filled () const { 1555 return filled_; 1556 } 1557 BOOST_UBLAS_INLINE index_data() const1558 const index_array_type &index_data () const { 1559 return index_data_; 1560 } 1561 BOOST_UBLAS_INLINE value_data() const1562 const value_array_type &value_data () const { 1563 return value_data_; 1564 } 1565 BOOST_UBLAS_INLINE set_filled(const typename index_array_type::size_type & sorted,const typename index_array_type::size_type & filled)1566 void set_filled (const typename index_array_type::size_type &sorted, const typename index_array_type::size_type &filled) { 1567 sorted_filled_ = sorted; 1568 filled_ = filled; 1569 storage_invariants (); 1570 } 1571 BOOST_UBLAS_INLINE index_data()1572 index_array_type &index_data () { 1573 return index_data_; 1574 } 1575 BOOST_UBLAS_INLINE value_data()1576 value_array_type &value_data () { 1577 return value_data_; 1578 } 1579 1580 // Resizing 1581 private: 1582 BOOST_UBLAS_INLINE restrict_capacity(size_type non_zeros) const1583 size_type restrict_capacity (size_type non_zeros) const { 1584 // minimum non_zeros 1585 non_zeros = (std::max) (non_zeros, size_type (1)); 1586 // ISSUE no maximum as coordinate may contain inserted duplicates 1587 return non_zeros; 1588 } 1589 public: 1590 BOOST_UBLAS_INLINE resize(size_type size,bool preserve=true)1591 void resize (size_type size, bool preserve = true) { 1592 if (preserve) 1593 sort (); // remove duplicate elements. 1594 size_ = size; 1595 capacity_ = restrict_capacity (capacity_); 1596 if (preserve) { 1597 index_data_. resize (capacity_, size_type ()); 1598 value_data_. resize (capacity_, value_type ()); 1599 filled_ = (std::min) (capacity_, filled_); 1600 while ((filled_ > 0) && (zero_based(index_data_[filled_ - 1]) >= size)) { 1601 --filled_; 1602 } 1603 } 1604 else { 1605 index_data_. resize (capacity_); 1606 value_data_. resize (capacity_); 1607 filled_ = 0; 1608 } 1609 sorted_filled_ = filled_; 1610 storage_invariants (); 1611 } 1612 // Reserving 1613 BOOST_UBLAS_INLINE reserve(size_type non_zeros,bool preserve=true)1614 void reserve (size_type non_zeros, bool preserve = true) { 1615 if (preserve) 1616 sort (); // remove duplicate elements. 1617 capacity_ = restrict_capacity (non_zeros); 1618 if (preserve) { 1619 index_data_. resize (capacity_, size_type ()); 1620 value_data_. resize (capacity_, value_type ()); 1621 filled_ = (std::min) (capacity_, filled_); 1622 } 1623 else { 1624 index_data_. resize (capacity_); 1625 value_data_. resize (capacity_); 1626 filled_ = 0; 1627 } 1628 sorted_filled_ = filled_; 1629 storage_invariants (); 1630 } 1631 1632 // Element support 1633 BOOST_UBLAS_INLINE find_element(size_type i)1634 pointer find_element (size_type i) { 1635 return const_cast<pointer> (const_cast<const self_type&>(*this).find_element (i)); 1636 } 1637 BOOST_UBLAS_INLINE find_element(size_type i) const1638 const_pointer find_element (size_type i) const { 1639 sort (); 1640 const_subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1641 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 1642 return 0; 1643 return &value_data_ [it - index_data_.begin ()]; 1644 } 1645 1646 // Element access 1647 BOOST_UBLAS_INLINE operator ()(size_type i) const1648 const_reference operator () (size_type i) const { 1649 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1650 sort (); 1651 const_subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1652 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 1653 return zero_; 1654 return value_data_ [it - index_data_.begin ()]; 1655 } 1656 BOOST_UBLAS_INLINE ref(size_type i)1657 true_reference ref (size_type i) { 1658 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1659 sort (); 1660 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1661 if (it == index_data_.begin () + filled_ || *it != k_based (i)) 1662 return insert_element (i, value_type/*zero*/()); 1663 else 1664 return value_data_ [it - index_data_.begin ()]; 1665 } 1666 BOOST_UBLAS_INLINE operator ()(size_type i)1667 reference operator () (size_type i) { 1668 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE 1669 return ref (i); 1670 #else 1671 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1672 return reference (*this, i); 1673 #endif 1674 } 1675 1676 BOOST_UBLAS_INLINE operator [](size_type i) const1677 const_reference operator [] (size_type i) const { 1678 return (*this) (i); 1679 } 1680 BOOST_UBLAS_INLINE operator [](size_type i)1681 reference operator [] (size_type i) { 1682 return (*this) (i); 1683 } 1684 1685 // Element assignment 1686 BOOST_UBLAS_INLINE append_element(size_type i,const_reference t)1687 void append_element (size_type i, const_reference t) { 1688 if (filled_ >= capacity_) 1689 reserve (2 * filled_, true); 1690 BOOST_UBLAS_CHECK (filled_ < capacity_, internal_logic ()); 1691 index_data_ [filled_] = k_based (i); 1692 value_data_ [filled_] = t; 1693 ++ filled_; 1694 sorted_ = false; 1695 storage_invariants (); 1696 } 1697 BOOST_UBLAS_INLINE insert_element(size_type i,const_reference t)1698 true_reference insert_element (size_type i, const_reference t) { 1699 BOOST_UBLAS_CHECK (!find_element (i), bad_index ()); // duplicate element 1700 append_element (i, t); 1701 return value_data_ [filled_ - 1]; 1702 } 1703 BOOST_UBLAS_INLINE erase_element(size_type i)1704 void erase_element (size_type i) { 1705 sort (); 1706 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1707 typename std::iterator_traits<subiterator_type>::difference_type n = it - index_data_.begin (); 1708 if (filled_ > typename index_array_type::size_type (n) && *it == k_based (i)) { 1709 std::copy (it + 1, index_data_.begin () + filled_, it); 1710 typename value_array_type::iterator itt (value_data_.begin () + n); 1711 std::copy (itt + 1, value_data_.begin () + filled_, itt); 1712 -- filled_; 1713 sorted_filled_ = filled_; 1714 } 1715 storage_invariants (); 1716 } 1717 1718 // Zeroing 1719 BOOST_UBLAS_INLINE clear()1720 void clear () { 1721 filled_ = 0; 1722 sorted_filled_ = filled_; 1723 sorted_ = true; 1724 storage_invariants (); 1725 } 1726 1727 // Assignment 1728 BOOST_UBLAS_INLINE operator =(const coordinate_vector & v)1729 coordinate_vector &operator = (const coordinate_vector &v) { 1730 if (this != &v) { 1731 size_ = v.size_; 1732 capacity_ = v.capacity_; 1733 filled_ = v.filled_; 1734 sorted_filled_ = v.sorted_filled_; 1735 sorted_ = v.sorted_; 1736 index_data_ = v.index_data_; 1737 value_data_ = v.value_data_; 1738 } 1739 storage_invariants (); 1740 return *this; 1741 } 1742 template<class C> // Container assignment without temporary 1743 BOOST_UBLAS_INLINE operator =(const vector_container<C> & v)1744 coordinate_vector &operator = (const vector_container<C> &v) { 1745 resize (v ().size (), false); 1746 assign (v); 1747 return *this; 1748 } 1749 BOOST_UBLAS_INLINE assign_temporary(coordinate_vector & v)1750 coordinate_vector &assign_temporary (coordinate_vector &v) { 1751 swap (v); 1752 return *this; 1753 } 1754 template<class AE> 1755 BOOST_UBLAS_INLINE operator =(const vector_expression<AE> & ae)1756 coordinate_vector &operator = (const vector_expression<AE> &ae) { 1757 self_type temporary (ae, capacity_); 1758 return assign_temporary (temporary); 1759 } 1760 template<class AE> 1761 BOOST_UBLAS_INLINE assign(const vector_expression<AE> & ae)1762 coordinate_vector &assign (const vector_expression<AE> &ae) { 1763 vector_assign<scalar_assign> (*this, ae); 1764 return *this; 1765 } 1766 1767 // Computed assignment 1768 template<class AE> 1769 BOOST_UBLAS_INLINE operator +=(const vector_expression<AE> & ae)1770 coordinate_vector &operator += (const vector_expression<AE> &ae) { 1771 self_type temporary (*this + ae, capacity_); 1772 return assign_temporary (temporary); 1773 } 1774 template<class C> // Container assignment without temporary 1775 BOOST_UBLAS_INLINE operator +=(const vector_container<C> & v)1776 coordinate_vector &operator += (const vector_container<C> &v) { 1777 plus_assign (v); 1778 return *this; 1779 } 1780 template<class AE> 1781 BOOST_UBLAS_INLINE plus_assign(const vector_expression<AE> & ae)1782 coordinate_vector &plus_assign (const vector_expression<AE> &ae) { 1783 vector_assign<scalar_plus_assign> (*this, ae); 1784 return *this; 1785 } 1786 template<class AE> 1787 BOOST_UBLAS_INLINE operator -=(const vector_expression<AE> & ae)1788 coordinate_vector &operator -= (const vector_expression<AE> &ae) { 1789 self_type temporary (*this - ae, capacity_); 1790 return assign_temporary (temporary); 1791 } 1792 template<class C> // Container assignment without temporary 1793 BOOST_UBLAS_INLINE operator -=(const vector_container<C> & v)1794 coordinate_vector &operator -= (const vector_container<C> &v) { 1795 minus_assign (v); 1796 return *this; 1797 } 1798 template<class AE> 1799 BOOST_UBLAS_INLINE minus_assign(const vector_expression<AE> & ae)1800 coordinate_vector &minus_assign (const vector_expression<AE> &ae) { 1801 vector_assign<scalar_minus_assign> (*this, ae); 1802 return *this; 1803 } 1804 template<class AT> 1805 BOOST_UBLAS_INLINE operator *=(const AT & at)1806 coordinate_vector &operator *= (const AT &at) { 1807 vector_assign_scalar<scalar_multiplies_assign> (*this, at); 1808 return *this; 1809 } 1810 template<class AT> 1811 BOOST_UBLAS_INLINE operator /=(const AT & at)1812 coordinate_vector &operator /= (const AT &at) { 1813 vector_assign_scalar<scalar_divides_assign> (*this, at); 1814 return *this; 1815 } 1816 1817 // Swapping 1818 BOOST_UBLAS_INLINE swap(coordinate_vector & v)1819 void swap (coordinate_vector &v) { 1820 if (this != &v) { 1821 std::swap (size_, v.size_); 1822 std::swap (capacity_, v.capacity_); 1823 std::swap (filled_, v.filled_); 1824 std::swap (sorted_filled_, v.sorted_filled_); 1825 std::swap (sorted_, v.sorted_); 1826 index_data_.swap (v.index_data_); 1827 value_data_.swap (v.value_data_); 1828 } 1829 storage_invariants (); 1830 } 1831 BOOST_UBLAS_INLINE swap(coordinate_vector & v1,coordinate_vector & v2)1832 friend void swap (coordinate_vector &v1, coordinate_vector &v2) { 1833 v1.swap (v2); 1834 } 1835 1836 // replacement if STL lower bound algorithm for use of inplace_merge lower_bound(size_type beg,size_type end,size_type target) const1837 size_type lower_bound (size_type beg, size_type end, size_type target) const { 1838 while (end > beg) { 1839 size_type mid = (beg + end) / 2; 1840 if (index_data_[mid] < index_data_[target]) { 1841 beg = mid + 1; 1842 } else { 1843 end = mid; 1844 } 1845 } 1846 return beg; 1847 } 1848 1849 // specialized replacement of STL inplace_merge to avoid compilation 1850 // problems with respect to the array_triple iterator inplace_merge(size_type beg,size_type mid,size_type end) const1851 void inplace_merge (size_type beg, size_type mid, size_type end) const { 1852 size_type len_lef = mid - beg; 1853 size_type len_rig = end - mid; 1854 1855 if (len_lef == 1 && len_rig == 1) { 1856 if (index_data_[mid] < index_data_[beg]) { 1857 std::swap(index_data_[beg], index_data_[mid]); 1858 std::swap(value_data_[beg], value_data_[mid]); 1859 } 1860 } else if (len_lef > 0 && len_rig > 0) { 1861 size_type lef_mid, rig_mid; 1862 if (len_lef >= len_rig) { 1863 lef_mid = (beg + mid) / 2; 1864 rig_mid = lower_bound(mid, end, lef_mid); 1865 } else { 1866 rig_mid = (mid + end) / 2; 1867 lef_mid = lower_bound(beg, mid, rig_mid); 1868 } 1869 std::rotate(&index_data_[0] + lef_mid, &index_data_[0] + mid, &index_data_[0] + rig_mid); 1870 std::rotate(&value_data_[0] + lef_mid, &value_data_[0] + mid, &value_data_[0] + rig_mid); 1871 1872 size_type new_mid = lef_mid + rig_mid - mid; 1873 inplace_merge(beg, lef_mid, new_mid); 1874 inplace_merge(new_mid, rig_mid, end); 1875 } 1876 } 1877 1878 // Sorting and summation of duplicates 1879 BOOST_UBLAS_INLINE sort() const1880 void sort () const { 1881 if (! sorted_ && filled_ > 0) { 1882 typedef index_pair_array<index_array_type, value_array_type> array_pair; 1883 array_pair ipa (filled_, index_data_, value_data_); 1884 #ifndef BOOST_UBLAS_COO_ALWAYS_DO_FULL_SORT 1885 const typename array_pair::iterator iunsorted = ipa.begin () + sorted_filled_; 1886 // sort new elements and merge 1887 std::sort (iunsorted, ipa.end ()); 1888 inplace_merge(0, sorted_filled_, filled_); 1889 #else 1890 const typename array_pair::iterator iunsorted = ipa.begin (); 1891 std::sort (iunsorted, ipa.end ()); 1892 #endif 1893 1894 // sum duplicates with += and remove 1895 size_type filled = 0; 1896 for (size_type i = 1; i < filled_; ++ i) { 1897 if (index_data_ [filled] != index_data_ [i]) { 1898 ++ filled; 1899 if (filled != i) { 1900 index_data_ [filled] = index_data_ [i]; 1901 value_data_ [filled] = value_data_ [i]; 1902 } 1903 } else { 1904 value_data_ [filled] += value_data_ [i]; 1905 } 1906 } 1907 filled_ = filled + 1; 1908 sorted_filled_ = filled_; 1909 sorted_ = true; 1910 storage_invariants (); 1911 } 1912 } 1913 1914 // Back element insertion and erasure 1915 BOOST_UBLAS_INLINE push_back(size_type i,const_reference t)1916 void push_back (size_type i, const_reference t) { 1917 // must maintain sort order 1918 BOOST_UBLAS_CHECK (sorted_ && (filled_ == 0 || index_data_ [filled_ - 1] < k_based (i)), external_logic ()); 1919 if (filled_ >= capacity_) 1920 reserve (2 * filled_, true); 1921 BOOST_UBLAS_CHECK (filled_ < capacity_, internal_logic ()); 1922 index_data_ [filled_] = k_based (i); 1923 value_data_ [filled_] = t; 1924 ++ filled_; 1925 sorted_filled_ = filled_; 1926 storage_invariants (); 1927 } 1928 BOOST_UBLAS_INLINE pop_back()1929 void pop_back () { 1930 // ISSUE invariants could be simpilfied if sorted required as precondition 1931 BOOST_UBLAS_CHECK (filled_ > 0, external_logic ()); 1932 -- filled_; 1933 sorted_filled_ = (std::min) (sorted_filled_, filled_); 1934 sorted_ = sorted_filled_ = filled_; 1935 storage_invariants (); 1936 } 1937 1938 // Iterator types 1939 private: 1940 // Use index array iterator 1941 typedef typename IA::const_iterator const_subiterator_type; 1942 typedef typename IA::iterator subiterator_type; 1943 1944 BOOST_UBLAS_INLINE at_element(size_type i)1945 true_reference at_element (size_type i) { 1946 BOOST_UBLAS_CHECK (i < size_, bad_index ()); 1947 sort (); 1948 subiterator_type it (detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1949 BOOST_UBLAS_CHECK (it != index_data_.begin () + filled_ && *it == k_based (i), bad_index ()); 1950 return value_data_ [it - index_data_.begin ()]; 1951 } 1952 1953 public: 1954 class const_iterator; 1955 class iterator; 1956 1957 // Element lookup 1958 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i) const1959 const_iterator find (size_type i) const { 1960 sort (); 1961 return const_iterator (*this, detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1962 } 1963 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. find(size_type i)1964 iterator find (size_type i) { 1965 sort (); 1966 return iterator (*this, detail::lower_bound (index_data_.begin (), index_data_.begin () + filled_, k_based (i), std::less<size_type> ())); 1967 } 1968 1969 1970 class const_iterator: 1971 public container_const_reference<coordinate_vector>, 1972 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 1973 const_iterator, value_type> { 1974 public: 1975 typedef typename coordinate_vector::value_type value_type; 1976 typedef typename coordinate_vector::difference_type difference_type; 1977 typedef typename coordinate_vector::const_reference reference; 1978 typedef const typename coordinate_vector::pointer pointer; 1979 1980 // Construction and destruction 1981 BOOST_UBLAS_INLINE const_iterator()1982 const_iterator (): 1983 container_const_reference<self_type> (), it_ () {} 1984 BOOST_UBLAS_INLINE const_iterator(const self_type & v,const const_subiterator_type & it)1985 const_iterator (const self_type &v, const const_subiterator_type &it): 1986 container_const_reference<self_type> (v), it_ (it) {} 1987 BOOST_UBLAS_INLINE const_iterator(const typename self_type::iterator & it)1988 const_iterator (const typename self_type::iterator &it): // ISSUE self_type:: stops VC8 using std::iterator here 1989 container_const_reference<self_type> (it ()), it_ (it.it_) {} 1990 1991 // Arithmetic 1992 BOOST_UBLAS_INLINE operator ++()1993 const_iterator &operator ++ () { 1994 ++ it_; 1995 return *this; 1996 } 1997 BOOST_UBLAS_INLINE operator --()1998 const_iterator &operator -- () { 1999 -- it_; 2000 return *this; 2001 } 2002 2003 // Dereference 2004 BOOST_UBLAS_INLINE operator *() const2005 const_reference operator * () const { 2006 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 2007 return (*this) ().value_data_ [it_ - (*this) ().index_data_.begin ()]; 2008 } 2009 2010 // Index 2011 BOOST_UBLAS_INLINE index() const2012 size_type index () const { 2013 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 2014 BOOST_UBLAS_CHECK ((*this) ().zero_based (*it_) < (*this) ().size (), bad_index ()); 2015 return (*this) ().zero_based (*it_); 2016 } 2017 2018 // Assignment 2019 BOOST_UBLAS_INLINE operator =(const const_iterator & it)2020 const_iterator &operator = (const const_iterator &it) { 2021 container_const_reference<self_type>::assign (&it ()); 2022 it_ = it.it_; 2023 return *this; 2024 } 2025 2026 // Comparison 2027 BOOST_UBLAS_INLINE operator ==(const const_iterator & it) const2028 bool operator == (const const_iterator &it) const { 2029 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 2030 return it_ == it.it_; 2031 } 2032 2033 private: 2034 const_subiterator_type it_; 2035 }; 2036 2037 BOOST_UBLAS_INLINE begin() const2038 const_iterator begin () const { 2039 return find (0); 2040 } 2041 BOOST_UBLAS_INLINE cbegin() const2042 const_iterator cbegin () const { 2043 return begin(); 2044 } 2045 BOOST_UBLAS_INLINE end() const2046 const_iterator end () const { 2047 return find (size_); 2048 } 2049 BOOST_UBLAS_INLINE cend() const2050 const_iterator cend () const { 2051 return end(); 2052 } 2053 2054 class iterator: 2055 public container_reference<coordinate_vector>, 2056 public bidirectional_iterator_base<sparse_bidirectional_iterator_tag, 2057 iterator, value_type> { 2058 public: 2059 typedef typename coordinate_vector::value_type value_type; 2060 typedef typename coordinate_vector::difference_type difference_type; 2061 typedef typename coordinate_vector::true_reference reference; 2062 typedef typename coordinate_vector::pointer pointer; 2063 2064 // Construction and destruction 2065 BOOST_UBLAS_INLINE iterator()2066 iterator (): 2067 container_reference<self_type> (), it_ () {} 2068 BOOST_UBLAS_INLINE iterator(self_type & v,const subiterator_type & it)2069 iterator (self_type &v, const subiterator_type &it): 2070 container_reference<self_type> (v), it_ (it) {} 2071 2072 // Arithmetic 2073 BOOST_UBLAS_INLINE operator ++()2074 iterator &operator ++ () { 2075 ++ it_; 2076 return *this; 2077 } 2078 BOOST_UBLAS_INLINE operator --()2079 iterator &operator -- () { 2080 -- it_; 2081 return *this; 2082 } 2083 2084 // Dereference 2085 BOOST_UBLAS_INLINE operator *() const2086 reference operator * () const { 2087 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ()); 2088 return (*this) ().value_data_ [it_ - (*this) ().index_data_.begin ()]; 2089 } 2090 2091 // Index 2092 BOOST_UBLAS_INLINE index() const2093 size_type index () const { 2094 BOOST_UBLAS_CHECK (*this != (*this) ().end (), bad_index ()); 2095 BOOST_UBLAS_CHECK ((*this) ().zero_based (*it_) < (*this) ().size (), bad_index ()); 2096 return (*this) ().zero_based (*it_); 2097 } 2098 2099 // Assignment 2100 BOOST_UBLAS_INLINE operator =(const iterator & it)2101 iterator &operator = (const iterator &it) { 2102 container_reference<self_type>::assign (&it ()); 2103 it_ = it.it_; 2104 return *this; 2105 } 2106 2107 // Comparison 2108 BOOST_UBLAS_INLINE operator ==(const iterator & it) const2109 bool operator == (const iterator &it) const { 2110 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ()); 2111 return it_ == it.it_; 2112 } 2113 2114 private: 2115 subiterator_type it_; 2116 2117 friend class const_iterator; 2118 }; 2119 2120 BOOST_UBLAS_INLINE begin()2121 iterator begin () { 2122 return find (0); 2123 } 2124 BOOST_UBLAS_INLINE end()2125 iterator end () { 2126 return find (size_); 2127 } 2128 2129 // Reverse iterator 2130 typedef reverse_iterator_base<const_iterator> const_reverse_iterator; 2131 typedef reverse_iterator_base<iterator> reverse_iterator; 2132 2133 BOOST_UBLAS_INLINE rbegin() const2134 const_reverse_iterator rbegin () const { 2135 return const_reverse_iterator (end ()); 2136 } 2137 BOOST_UBLAS_INLINE crbegin() const2138 const_reverse_iterator crbegin () const { 2139 return rbegin (); 2140 } 2141 BOOST_UBLAS_INLINE rend() const2142 const_reverse_iterator rend () const { 2143 return const_reverse_iterator (begin ()); 2144 } 2145 BOOST_UBLAS_INLINE crend() const2146 const_reverse_iterator crend () const { 2147 return rend (); 2148 } 2149 BOOST_UBLAS_INLINE rbegin()2150 reverse_iterator rbegin () { 2151 return reverse_iterator (end ()); 2152 } 2153 BOOST_UBLAS_INLINE rend()2154 reverse_iterator rend () { 2155 return reverse_iterator (begin ()); 2156 } 2157 2158 // Serialization 2159 template<class Archive> serialize(Archive & ar,const unsigned int)2160 void serialize(Archive & ar, const unsigned int /* file_version */){ 2161 serialization::collection_size_type s (size_); 2162 ar & serialization::make_nvp("size",s); 2163 if (Archive::is_loading::value) { 2164 size_ = s; 2165 } 2166 // ISSUE: filled may be much less than capacity 2167 // ISSUE: index_data_ and value_data_ are undefined between filled and capacity (trouble with 'nan'-values) 2168 ar & serialization::make_nvp("capacity", capacity_); 2169 ar & serialization::make_nvp("filled", filled_); 2170 ar & serialization::make_nvp("sorted_filled", sorted_filled_); 2171 ar & serialization::make_nvp("sorted", sorted_); 2172 ar & serialization::make_nvp("index_data", index_data_); 2173 ar & serialization::make_nvp("value_data", value_data_); 2174 storage_invariants(); 2175 } 2176 2177 private: storage_invariants() const2178 void storage_invariants () const 2179 { 2180 BOOST_UBLAS_CHECK (capacity_ == index_data_.size (), internal_logic ()); 2181 BOOST_UBLAS_CHECK (capacity_ == value_data_.size (), internal_logic ()); 2182 BOOST_UBLAS_CHECK (filled_ <= capacity_, internal_logic ()); 2183 BOOST_UBLAS_CHECK (sorted_filled_ <= filled_, internal_logic ()); 2184 BOOST_UBLAS_CHECK (sorted_ == (sorted_filled_ == filled_), internal_logic ()); 2185 BOOST_UBLAS_CHECK ((0 == filled_) || (zero_based(index_data_[filled_ - 1]) < size_), internal_logic ()); 2186 } 2187 2188 size_type size_; 2189 size_type capacity_; 2190 mutable typename index_array_type::size_type filled_; 2191 mutable typename index_array_type::size_type sorted_filled_; 2192 mutable bool sorted_; 2193 mutable index_array_type index_data_; 2194 mutable value_array_type value_data_; 2195 static const value_type zero_; 2196 2197 BOOST_UBLAS_INLINE zero_based(size_type k_based_index)2198 static size_type zero_based (size_type k_based_index) { 2199 return k_based_index - IB; 2200 } 2201 BOOST_UBLAS_INLINE k_based(size_type zero_based_index)2202 static size_type k_based (size_type zero_based_index) { 2203 return zero_based_index + IB; 2204 } 2205 2206 friend class iterator; 2207 friend class const_iterator; 2208 }; 2209 2210 template<class T, std::size_t IB, class IA, class TA> 2211 const typename coordinate_vector<T, IB, IA, TA>::value_type coordinate_vector<T, IB, IA, TA>::zero_ = value_type/*zero*/(); 2212 2213 }}} 2214 2215 #endif 2216