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