boost/numeric/ublas/vector_of_vector.hpp

   1 //
   2 //  Copyright (c) 2003
   3 //  Gunter Winkler, Joerg Walter
   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_OF_VECTOR_
  14 #define _BOOST_UBLAS_VECTOR_OF_VECTOR_
  15
  16 #include <boost/type_traits.hpp>
  17
  18 #include <boost/numeric/ublas/storage_sparse.hpp>
  19 #include <boost/numeric/ublas/matrix_sparse.hpp>
  20
  21 // Iterators based on ideas of Jeremy Siek
  22
  23 namespace boost { namespace numeric { namespace ublas {
  24
  25     // uBLAS sparse vector based sparse matrix class
  26     // FIXME outer vector can be sparse type but it is completely filled
  27     template<class T, class L, class A>
  28     class generalized_vector_of_vector:
  29         public matrix_container<generalized_vector_of_vector<T, L, A> > {
  30
  31         typedef T &true_reference;
  32         typedef T *pointer;
  33         typedef const T *const_pointer;
  34         typedef L layout_type;
  35         typedef generalized_vector_of_vector<T, L, A> self_type;
  36     public:
  37 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
  38         using matrix_container<self_type>::operator ();
  39 #endif
  40         typedef typename A::size_type size_type;
  41         typedef typename A::difference_type difference_type;
  42         typedef T value_type;
  43         typedef const T &const_reference;
  44 #ifndef BOOST_UBLAS_STRICT_VECTOR_SPARSE
  45         typedef T &reference;
  46 #else
  47         typedef sparse_matrix_element<self_type> reference;
  48 #endif
  49         typedef A array_type;
  50         typedef const matrix_reference<const self_type> const_closure_type;
  51         typedef matrix_reference<self_type> closure_type;
  52         typedef typename A::value_type vector_data_value_type;
  53         typedef vector_data_value_type vector_temporary_type;
  54         typedef self_type matrix_temporary_type;
  55         typedef sparse_tag storage_category;
  56         typedef typename L::orientation_category orientation_category;
  57
  58         // Construction and destruction
  59         BOOST_UBLAS_INLINE
  60         generalized_vector_of_vector ():
  61             matrix_container<self_type> (),
  62             size1_ (0), size2_ (0), data_ (1) {
  63             const size_type sizeM = layout_type::size_M (size1_, size2_);
  64              // create size1+1 empty vector elements
  65             data_.insert_element (sizeM, vector_data_value_type ());
  66             storage_invariants ();
  67         }
  68         BOOST_UBLAS_INLINE
  69         generalized_vector_of_vector (size_type size1, size_type size2, size_type non_zeros = 0):
  70             matrix_container<self_type> (),
  71             size1_ (size1), size2_ (size2), data_ (layout_type::size_M (size1_, size2_) + 1) {
  72             const size_type sizeM = layout_type::size_M (size1_, size2_);
  73             const size_type sizem = layout_type::size_m (size1_, size2_);
  74             for (size_type i = 0; i < sizeM; ++ i) // create size1 vector elements
  75                 data_.insert_element (i, vector_data_value_type ()) .resize (sizem, false);
  76             data_.insert_element (sizeM, vector_data_value_type ());
  77             storage_invariants ();
  78         }
  79         BOOST_UBLAS_INLINE
  80         generalized_vector_of_vector (const generalized_vector_of_vector &m):
  81             matrix_container<self_type> (),
  82             size1_ (m.size1_), size2_ (m.size2_), data_ (m.data_) {
  83             storage_invariants ();
  84         }
  85         template<class AE>
  86         BOOST_UBLAS_INLINE
  87         generalized_vector_of_vector (const matrix_expression<AE> &ae, size_type non_zeros = 0):
  88             matrix_container<self_type> (),
  89             size1_ (ae ().size1 ()), size2_ (ae ().size2 ()), data_ (layout_type::size_M (size1_, size2_) + 1) {
  90             const size_type sizeM = layout_type::size_M (size1_, size2_);
  91             const size_type sizem = layout_type::size_m (size1_, size2_);
  92             for (size_type i = 0; i < sizeM; ++ i) // create size1 vector elements
  93                 data_.insert_element (i, vector_data_value_type ()) .resize (sizem, false);
  94             data_.insert_element (sizeM, vector_data_value_type ());
  95             storage_invariants ();
  96             matrix_assign<scalar_assign> (*this, ae);
  97         }
  98
  99         // Accessors
 100         BOOST_UBLAS_INLINE
 101         size_type size1 () const {
 102             return size1_;
 103         }
 104         BOOST_UBLAS_INLINE
 105         size_type size2 () const {
 106             return size2_;
 107         }
 108         BOOST_UBLAS_INLINE
 109         size_type nnz_capacity () const {
 110             size_type non_zeros = 0;
 111             for (const_vectoriterator_type itv = data_.begin (); itv != data_.end (); ++ itv)
 112                 non_zeros += (*itv).nnz_capacity ();
 113             return non_zeros;
 114         }
 115         BOOST_UBLAS_INLINE
 116         size_type nnz () const {
 117             size_type non_zeros = 0;
 118             for (const_vectoriterator_type itv = data_.begin (); itv != data_.end (); ++ itv)
 119                 non_zeros += (*itv).nnz ();
 120             return non_zeros;
 121         }
 122
 123         // Storage accessors
 124         BOOST_UBLAS_INLINE
 125         const array_type &data () const {
 126             return data_;
 127         }
 128         BOOST_UBLAS_INLINE
 129         array_type &data () {
 130             return data_;
 131         }
 132
 133         // Resizing
 134         BOOST_UBLAS_INLINE
 135         void resize (size_type size1, size_type size2, bool preserve = true) {
 136             const size_type oldM = layout_type::size_M (size1_, size2_);
 137             size1_ = size1;
 138             size2_ = size2;
 139             const size_type sizeM = layout_type::size_M (size1_, size2_);
 140             const size_type sizem = layout_type::size_m (size1_, size2_);
 141             data ().resize (sizeM + 1, preserve);
 142             if (preserve) {
 143                 for (size_type i = 0; (i <= oldM) && (i < sizeM); ++ i)
 144                     ref (data () [i]).resize (sizem, preserve);
 145                 for (size_type i = oldM+1; i < sizeM; ++ i) // create new vector elements
 146                     data_.insert_element (i, vector_data_value_type ()) .resize (sizem, false);
 147                 if (sizeM > oldM) {
 148                     data_.insert_element (sizeM, vector_data_value_type ());
 149                 } else {
 150                     ref (data () [sizeM]).resize (0, false);
 151                 }
 152             } else {
 153                 for (size_type i = 0; i < sizeM; ++ i)
 154                     data_.insert_element (i, vector_data_value_type ()) .resize (sizem, false);
 155                 data_.insert_element (sizeM, vector_data_value_type ());
 156             }
 157             storage_invariants ();
 158         }
 159
 160         // Element support
 161         BOOST_UBLAS_INLINE
 162         pointer find_element (size_type i, size_type j) {
 163             return const_cast<pointer> (const_cast<const self_type&>(*this).find_element (i, j));
 164         }
 165         BOOST_UBLAS_INLINE
 166         const_pointer find_element (size_type i, size_type j) const {
 167             const size_type elementM = layout_type::index_M (i, j);
 168             const size_type elementm = layout_type::index_m (i, j);
 169             // optimise: check the storage_type and index directly if element always exists
 170             if (boost::is_convertible<typename array_type::storage_category, packed_tag>::value) {
 171                 return & (data () [elementM] [elementm]);
 172             }
 173             else {
 174                 const typename array_type::value_type *pv = data ().find_element (elementM);
 175                 if (!pv)
 176                     return 0;
 177                 return pv->find_element (elementm);
 178             }
 179         }
 180
 181         // Element access
 182         BOOST_UBLAS_INLINE
 183         const_reference operator () (size_type i, size_type j) const {
 184             const_pointer p = find_element (i, j);
 185             // optimise: check the storage_type and index directly if element always exists
 186             if (boost::is_convertible<typename array_type::storage_category, packed_tag>::value) {
 187                 BOOST_UBLAS_CHECK (p, internal_logic () );
 188                 return *p;
 189             }
 190             else {
 191                 if (p)
 192                     return *p;
 193                 else
 194                     return zero_;
 195             }
 196         }
 197         BOOST_UBLAS_INLINE
 198         reference operator () (size_type i, size_type j) {
 199 #ifndef BOOST_UBLAS_STRICT_MATRIX_SPARSE
 200             return at_element (i, j);
 201 #else
 202             return reference (*this, i, j);
 203 #endif
 204         }
 205
 206         // Assignment
 207         BOOST_UBLAS_INLINE
 208         generalized_vector_of_vector &operator = (const generalized_vector_of_vector &m) {
 209             if (this != &m) {
 210                 size1_ = m.size1_;
 211                 size2_ = m.size2_;
 212                 data () = m.data ();
 213             }
 214             storage_invariants ();
 215             return *this;
 216         }
 217         BOOST_UBLAS_INLINE
 218         generalized_vector_of_vector &assign_temporary (generalized_vector_of_vector &m) {
 219             swap (m);
 220             return *this;
 221         }
 222         template<class AE>
 223         BOOST_UBLAS_INLINE
 224         generalized_vector_of_vector &operator = (const matrix_expression<AE> &ae) {
 225             self_type temporary (ae);
 226             return assign_temporary (temporary);
 227         }
 228         template<class AE>
 229         BOOST_UBLAS_INLINE
 230         generalized_vector_of_vector &assign (const matrix_expression<AE> &ae) {
 231             matrix_assign<scalar_assign> (*this, ae);
 232             return *this;
 233         }
 234         template<class AE>
 235         BOOST_UBLAS_INLINE
 236         generalized_vector_of_vector& operator += (const matrix_expression<AE> &ae) {
 237             self_type temporary (*this + ae);
 238             return assign_temporary (temporary);
 239         }
 240         template<class AE>
 241         BOOST_UBLAS_INLINE
 242         generalized_vector_of_vector &plus_assign (const matrix_expression<AE> &ae) {
 243             matrix_assign<scalar_plus_assign> (*this, ae);
 244             return *this;
 245         }
 246         template<class AE>
 247         BOOST_UBLAS_INLINE
 248         generalized_vector_of_vector& operator -= (const matrix_expression<AE> &ae) {
 249             self_type temporary (*this - ae);
 250             return assign_temporary (temporary);
 251         }
 252         template<class AE>
 253         BOOST_UBLAS_INLINE
 254         generalized_vector_of_vector &minus_assign (const matrix_expression<AE> &ae) {
 255             matrix_assign<scalar_minus_assign> (*this, ae);
 256             return *this;
 257         }
 258         template<class AT>
 259         BOOST_UBLAS_INLINE
 260         generalized_vector_of_vector& operator *= (const AT &at) {
 261             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
 262             return *this;
 263         }
 264         template<class AT>
 265         BOOST_UBLAS_INLINE
 266         generalized_vector_of_vector& operator /= (const AT &at) {
 267             matrix_assign_scalar<scalar_divides_assign> (*this, at);
 268             return *this;
 269         }
 270
 271         // Swapping
 272         BOOST_UBLAS_INLINE
 273         void swap (generalized_vector_of_vector &m) {
 274             if (this != &m) {
 275                 std::swap (size1_, m.size1_);
 276                 std::swap (size2_, m.size2_);
 277                 data ().swap (m.data ());
 278             }
 279             storage_invariants ();
 280         }
 281         BOOST_UBLAS_INLINE
 282         friend void swap (generalized_vector_of_vector &m1, generalized_vector_of_vector &m2) {
 283             m1.swap (m2);
 284         }
 285
 286         // Sorting
 287         void sort () {
 288             vectoriterator_type itv (data ().begin ());
 289             vectoriterator_type itv_end (data ().end ());
 290             while (itv != itv_end) {
 291                 (*itv).sort ();
 292                 ++ itv;
 293             }
 294         }
 295
 296         // Element insertion and erasure
 297         BOOST_UBLAS_INLINE
 298         true_reference insert_element (size_type i, size_type j, const_reference t) {
 299             const size_type elementM = layout_type::index_M (i, j);
 300             const size_type elementm = layout_type::index_m (i, j);
 301             vector_data_value_type& vd (ref (data () [elementM]));
 302             storage_invariants ();
 303             return vd.insert_element (elementm, t);
 304         }
 305         BOOST_UBLAS_INLINE
 306         void append_element (size_type i, size_type j, const_reference t) {
 307             const size_type elementM = layout_type::index_M (i, j);
 308             const size_type elementm = layout_type::index_m (i, j);
 309             vector_data_value_type& vd (ref (data () [elementM]));
 310             storage_invariants ();
 311             return vd.append_element (elementm, t);
 312         }
 313         BOOST_UBLAS_INLINE
 314         void erase_element (size_type i, size_type j) {
 315             vectoriterator_type itv (data ().find (layout_type::index_M (i, j)));
 316             if (itv == data ().end ())
 317                 return;
 318             (*itv).erase_element (layout_type::index_m (i, j));
 319             storage_invariants ();
 320         }
 321         BOOST_UBLAS_INLINE
 322         void clear () {
 323             const size_type sizeM = layout_type::size_M (size1_, size2_);
 324             // FIXME should clear data () if this is done via value_type/*zero*/() then it is not size preserving
 325             for (size_type i = 0; i < sizeM; ++ i)
 326                 ref (data () [i]).clear ();
 327             storage_invariants ();
 328         }
 329
 330         // Iterator types
 331     private:
 332         // Use vector iterator
 333         typedef typename A::const_iterator const_vectoriterator_type;
 334         typedef typename A::iterator vectoriterator_type;
 335         typedef typename A::value_type::const_iterator const_subiterator_type;
 336         typedef typename A::value_type::iterator subiterator_type;
 337
 338         BOOST_UBLAS_INLINE
 339         true_reference at_element (size_type i, size_type j) {
 340             return ref (ref (data () [layout_type::index_M (i, j)]) [layout_type::index_m (i, j)]);
 341         }
 342
 343     public:
 344         class const_iterator1;
 345         class iterator1;
 346         class const_iterator2;
 347         class iterator2;
 348         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
 349         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
 350         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
 351         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
 352
 353         // Element lookup
 354         // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
 355         const_iterator1 find1 (int rank, size_type i, size_type j, int direction = 1) const {
 356             for (;;) {
 357                 const_vectoriterator_type itv (data ().find (layout_type::index_M (i, j)));
 358                 const_vectoriterator_type itv_end (data ().end ());
 359                 if (itv == itv_end)
 360                     return const_iterator1 (*this, rank, i, j, itv_end, (*(-- itv)).end ());
 361
 362                 const_subiterator_type it ((*itv).find (layout_type::index_m (i, j)));
 363                 const_subiterator_type it_end ((*itv).end ());
 364                 if (rank == 0)
 365                     return const_iterator1 (*this, rank, i, j, itv, it);
 366                 if (it != it_end && it.index () == layout_type::index_m (i, j))
 367                     return const_iterator1 (*this, rank, i, j, itv, it);
 368                 if (direction > 0) {
 369                     if (layout_type::fast_i ()) {
 370                         if (it == it_end)
 371                             return const_iterator1 (*this, rank, i, j, itv, it);
 372                         i = it.index ();
 373                     } else {
 374                         if (i >= size1_)
 375                             return const_iterator1 (*this, rank, i, j, itv, it);
 376                         ++ i;
 377                     }
 378                 } else /* if (direction < 0)  */ {
 379                     if (layout_type::fast_i ()) {
 380                         if (it == (*itv).begin ())
 381                             return const_iterator1 (*this, rank, i, j, itv, it);
 382                         --it;
 383                         i = it.index ();
 384                     } else {
 385                         if (i == 0)
 386                             return const_iterator1 (*this, rank, i, j, itv, it);
 387                         -- i;
 388                     }
 389                 }
 390             }
 391         }
 392         // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
 393         iterator1 find1 (int rank, size_type i, size_type j, int direction = 1) {
 394             for (;;) {
 395                 vectoriterator_type itv (data ().find (layout_type::index_M (i, j)));
 396                 vectoriterator_type itv_end (data ().end ());
 397                 if (itv == itv_end)
 398                     return iterator1 (*this, rank, i, j, itv_end, (*(-- itv)).end ());
 399
 400                 subiterator_type it ((*itv).find (layout_type::index_m (i, j)));
 401                 subiterator_type it_end ((*itv).end ());
 402                 if (rank == 0)
 403                     return iterator1 (*this, rank, i, j, itv, it);
 404                 if (it != it_end && it.index () == layout_type::index_m (i, j))
 405                     return iterator1 (*this, rank, i, j, itv, it);
 406                 if (direction > 0) {
 407                     if (layout_type::fast_i ()) {
 408                         if (it == it_end)
 409                             return iterator1 (*this, rank, i, j, itv, it);
 410                         i = it.index ();
 411                     } else {
 412                         if (i >= size1_)
 413                             return iterator1 (*this, rank, i, j, itv, it);
 414                         ++ i;
 415                     }
 416                 } else /* if (direction < 0)  */ {
 417                     if (layout_type::fast_i ()) {
 418                         if (it == (*itv).begin ())
 419                             return iterator1 (*this, rank, i, j, itv, it);
 420                         --it;
 421                         i = it.index ();
 422                     } else {
 423                         if (i == 0)
 424                             return iterator1 (*this, rank, i, j, itv, it);
 425                         -- i;
 426                     }
 427                 }
 428             }
 429         }
 430         // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
 431         const_iterator2 find2 (int rank, size_type i, size_type j, int direction = 1) const {
 432             for (;;) {
 433                 const_vectoriterator_type itv (data ().find (layout_type::index_M (i, j)));
 434                 const_vectoriterator_type itv_end (data ().end ());
 435                 if (itv == itv_end)
 436                     return const_iterator2 (*this, rank, i, j, itv_end, (*(-- itv)).end ());
 437
 438                 const_subiterator_type it ((*itv).find (layout_type::index_m (i, j)));
 439                 const_subiterator_type it_end ((*itv).end ());
 440                 if (rank == 0)
 441                     return const_iterator2 (*this, rank, i, j, itv, it);
 442                 if (it != it_end && it.index () == layout_type::index_m (i, j))
 443                     return const_iterator2 (*this, rank, i, j, itv, it);
 444                 if (direction > 0) {
 445                     if (layout_type::fast_j ()) {
 446                         if (it == it_end)
 447                             return const_iterator2 (*this, rank, i, j, itv, it);
 448                         j = it.index ();
 449                     } else {
 450                         if (j >= size2_)
 451                             return const_iterator2 (*this, rank, i, j, itv, it);
 452                         ++ j;
 453                     }
 454                 } else /* if (direction < 0)  */ {
 455                     if (layout_type::fast_j ()) {
 456                         if (it == (*itv).begin ())
 457                             return const_iterator2 (*this, rank, i, j, itv, it);
 458                         --it;
 459                         j = it.index ();
 460                     } else {
 461                         if (j == 0)
 462                             return const_iterator2 (*this, rank, i, j, itv, it);
 463                         -- j;
 464                     }
 465                 }
 466             }
 467         }
 468         // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
 469         iterator2 find2 (int rank, size_type i, size_type j, int direction = 1) {
 470             for (;;) {
 471                 vectoriterator_type itv (data ().find (layout_type::index_M (i, j)));
 472                 vectoriterator_type itv_end (data ().end ());
 473                 if (itv == itv_end)
 474                     return iterator2 (*this, rank, i, j, itv_end, (*(-- itv)).end ());
 475
 476                 subiterator_type it ((*itv).find (layout_type::index_m (i, j)));
 477                 subiterator_type it_end ((*itv).end ());
 478                 if (rank == 0)
 479                     return iterator2 (*this, rank, i, j, itv, it);
 480                 if (it != it_end && it.index () == layout_type::index_m (i, j))
 481                     return iterator2 (*this, rank, i, j, itv, it);
 482                 if (direction > 0) {
 483                     if (layout_type::fast_j ()) {
 484                         if (it == it_end)
 485                             return iterator2 (*this, rank, i, j, itv, it);
 486                         j = it.index ();
 487                     } else {
 488                         if (j >= size2_)
 489                             return iterator2 (*this, rank, i, j, itv, it);
 490                         ++ j;
 491                     }
 492                 } else /* if (direction < 0)  */ {
 493                     if (layout_type::fast_j ()) {
 494                         if (it == (*itv).begin ())
 495                             return iterator2 (*this, rank, i, j, itv, it);
 496                         --it;
 497                         j = it.index ();
 498                     } else {
 499                         if (j == 0)
 500                             return iterator2 (*this, rank, i, j, itv, it);
 501                         -- j;
 502                     }
 503                 }
 504             }
 505         }
 506
 507
 508         class const_iterator1:
 509             public container_const_reference<generalized_vector_of_vector>,
 510             public bidirectional_iterator_base<sparse_bidirectional_iterator_tag,
 511                                                const_iterator1, value_type> {
 512         public:
 513             typedef typename generalized_vector_of_vector::difference_type difference_type;
 514             typedef typename generalized_vector_of_vector::value_type value_type;
 515             typedef typename generalized_vector_of_vector::const_reference reference;
 516             typedef const typename generalized_vector_of_vector::pointer pointer;
 517
 518             typedef const_iterator2 dual_iterator_type;
 519             typedef const_reverse_iterator2 dual_reverse_iterator_type;
 520
 521             // Construction and destruction
 522             BOOST_UBLAS_INLINE
 523             const_iterator1 ():
 524                 container_const_reference<self_type> (), rank_ (), i_ (), j_ (), itv_ (), it_ () {}
 525             BOOST_UBLAS_INLINE
 526             const_iterator1 (const self_type &m, int rank, size_type i, size_type j, const const_vectoriterator_type &itv, const const_subiterator_type &it):
 527                 container_const_reference<self_type> (m), rank_ (rank), i_ (i), j_ (j), itv_ (itv), it_ (it) {}
 528             BOOST_UBLAS_INLINE
 529             const_iterator1 (const iterator1 &it):
 530                 container_const_reference<self_type> (it ()), rank_ (it.rank_), i_ (it.i_), j_ (it.j_), itv_ (it.itv_), it_ (it.it_) {}
 531
 532             // Arithmetic
 533             BOOST_UBLAS_INLINE
 534             const_iterator1 &operator ++ () {
 535                 if (rank_ == 1 && layout_type::fast_i ())
 536                     ++ it_;
 537                 else {
 538                     const self_type &m = (*this) ();
 539                     i_ = index1 () + 1;
 540                     if (rank_ == 1 && ++ itv_ == m.end1 ().itv_)
 541                         *this = m.find1 (rank_, i_, j_, 1);
 542                     else if (rank_ == 1) {
 543                         it_ = (*itv_).begin ();
 544                         if (it_ == (*itv_).end () || index2 () != j_)
 545                             *this = m.find1 (rank_, i_, j_, 1);
 546                     }
 547                 }
 548                 return *this;
 549             }
 550             BOOST_UBLAS_INLINE
 551             const_iterator1 &operator -- () {
 552                 if (rank_ == 1 && layout_type::fast_i ())
 553                     -- it_;
 554                 else {
 555                     const self_type &m = (*this) ();
 556                     i_ = index1 () - 1;
 557                     if (rank_ == 1 && -- itv_ == m.end1 ().itv_)
 558                         *this = m.find1 (rank_, i_, j_, -1);
 559                     else if (rank_ == 1) {
 560                         it_ = (*itv_).begin ();
 561                         if (it_ == (*itv_).end () || index2 () != j_)
 562                             *this = m.find1 (rank_, i_, j_, -1);
 563                     }
 564                 }
 565                 return *this;
 566             }
 567
 568             // Dereference
 569             BOOST_UBLAS_INLINE
 570             const_reference operator * () const {
 571                 BOOST_UBLAS_CHECK (index1 () < (*this) ().size1 (), bad_index ());
 572                 BOOST_UBLAS_CHECK (index2 () < (*this) ().size2 (), bad_index ());
 573                 if (rank_ == 1) {
 574                     return *it_;
 575                 } else {
 576                     return (*this) () (i_, j_);
 577                 }
 578             }
 579
 580 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
 581             BOOST_UBLAS_INLINE
 582 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 583             typename self_type::
 584 #endif
 585             const_iterator2 begin () const {
 586                 const self_type &m = (*this) ();
 587                 return m.find2 (1, index1 (), 0);
 588             }
 589             BOOST_UBLAS_INLINE
 590 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 591             typename self_type::
 592 #endif
 593             const_iterator2 end () const {
 594                 const self_type &m = (*this) ();
 595                 return m.find2 (1, index1 (), m.size2 ());
 596             }
 597             BOOST_UBLAS_INLINE
 598 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 599             typename self_type::
 600 #endif
 601             const_reverse_iterator2 rbegin () const {
 602                 return const_reverse_iterator2 (end ());
 603             }
 604             BOOST_UBLAS_INLINE
 605 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 606             typename self_type::
 607 #endif
 608             const_reverse_iterator2 rend () const {
 609                 return const_reverse_iterator2 (begin ());
 610             }
 611 #endif
 612
 613             // Indices
 614             BOOST_UBLAS_INLINE
 615             size_type index1 () const {
 616                 BOOST_UBLAS_CHECK (*this != (*this) ().find1 (0, (*this) ().size1 (), j_), bad_index ());
 617                 if (rank_ == 1) {
 618                     BOOST_UBLAS_CHECK (layout_type::index_M (itv_.index (), it_.index ()) < (*this) ().size1 (), bad_index ());
 619                     return layout_type::index_M (itv_.index (), it_.index ());
 620                 } else {
 621                     return i_;
 622                 }
 623             }
 624             BOOST_UBLAS_INLINE
 625             size_type index2 () const {
 626                 BOOST_UBLAS_CHECK (*this != (*this) ().find1 (0, (*this) ().size1 (), j_), bad_index ());
 627                 if (rank_ == 1) {
 628                     BOOST_UBLAS_CHECK (layout_type::index_m (itv_.index (), it_.index ()) < (*this) ().size2 (), bad_index ());
 629                     return layout_type::index_m (itv_.index (), it_.index ());
 630                 } else {
 631                     return j_;
 632                 }
 633             }
 634
 635             // Assignment
 636             BOOST_UBLAS_INLINE
 637             const_iterator1 &operator = (const const_iterator1 &it) {
 638                 container_const_reference<self_type>::assign (&it ());
 639                 rank_ = it.rank_;
 640                 i_ = it.i_;
 641                 j_ = it.j_;
 642                 itv_ = it.itv_;
 643                 it_ = it.it_;
 644                 return *this;
 645             }
 646
 647             // Comparison
 648             BOOST_UBLAS_INLINE
 649             bool operator == (const const_iterator1 &it) const {
 650                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
 651                 // BOOST_UBLAS_CHECK (rank_ == it.rank_, internal_logic ());
 652                 if (rank_ == 1 || it.rank_ == 1) {
 653                     return it_ == it.it_;
 654                 } else {
 655                     return i_ == it.i_ && j_ == it.j_;
 656                 }
 657             }
 658
 659         private:
 660             int rank_;
 661             size_type i_;
 662             size_type j_;
 663             const_vectoriterator_type itv_;
 664             const_subiterator_type it_;
 665         };
 666
 667         BOOST_UBLAS_INLINE
 668         const_iterator1 begin1 () const {
 669             return find1 (0, 0, 0);
 670         }
 671         BOOST_UBLAS_INLINE
 672         const_iterator1 end1 () const {
 673             return find1 (0, size1_, 0);
 674         }
 675
 676         class iterator1:
 677             public container_reference<generalized_vector_of_vector>,
 678             public bidirectional_iterator_base<sparse_bidirectional_iterator_tag,
 679                                                iterator1, value_type> {
 680         public:
 681             typedef typename generalized_vector_of_vector::difference_type difference_type;
 682             typedef typename generalized_vector_of_vector::value_type value_type;
 683             typedef typename generalized_vector_of_vector::true_reference reference;
 684             typedef typename generalized_vector_of_vector::pointer pointer;
 685
 686             typedef iterator2 dual_iterator_type;
 687             typedef reverse_iterator2 dual_reverse_iterator_type;
 688
 689             // Construction and destruction
 690             BOOST_UBLAS_INLINE
 691             iterator1 ():
 692                 container_reference<self_type> (), rank_ (), i_ (), j_ (), itv_ (), it_ () {}
 693             BOOST_UBLAS_INLINE
 694             iterator1 (self_type &m, int rank, size_type i, size_type j, const vectoriterator_type &itv, const subiterator_type &it):
 695                 container_reference<self_type> (m), rank_ (rank), i_ (i), j_ (j), itv_ (itv), it_ (it) {}
 696
 697             // Arithmetic
 698             BOOST_UBLAS_INLINE
 699             iterator1 &operator ++ () {
 700                 if (rank_ == 1 && layout_type::fast_i ())
 701                     ++ it_;
 702                 else {
 703                     self_type &m = (*this) ();
 704                     i_ = index1 () + 1;
 705                     if (rank_ == 1 && ++ itv_ == m.end1 ().itv_)
 706                         *this = m.find1 (rank_, i_, j_, 1);
 707                     else if (rank_ == 1) {
 708                         it_ = (*itv_).begin ();
 709                         if (it_ == (*itv_).end () || index2 () != j_)
 710                             *this = m.find1 (rank_, i_, j_, 1);
 711                     }
 712                 }
 713                 return *this;
 714             }
 715             BOOST_UBLAS_INLINE
 716             iterator1 &operator -- () {
 717                 if (rank_ == 1 && layout_type::fast_i ())
 718                     -- it_;
 719                 else {
 720                     self_type &m = (*this) ();
 721                     i_ = index1 () - 1;
 722                     if (rank_ == 1 && -- itv_ == m.end1 ().itv_)
 723                         *this = m.find1 (rank_, i_, j_, -1);
 724                     else if (rank_ == 1) {
 725                         it_ = (*itv_).begin ();
 726                         if (it_ == (*itv_).end () || index2 () != j_)
 727                             *this = m.find1 (rank_, i_, j_, -1);
 728                     }
 729                 }
 730                 return *this;
 731             }
 732
 733             // Dereference
 734             BOOST_UBLAS_INLINE
 735             true_reference operator * () const {
 736                 BOOST_UBLAS_CHECK (index1 () < (*this) ().size1 (), bad_index ());
 737                 BOOST_UBLAS_CHECK (index2 () < (*this) ().size2 (), bad_index ());
 738                 if (rank_ == 1) {
 739                     return *it_;
 740                 } else {
 741                     return (*this) ().at_element (i_, j_);
 742                 }
 743             }
 744
 745 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
 746             BOOST_UBLAS_INLINE
 747 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 748             typename self_type::
 749 #endif
 750             iterator2 begin () const {
 751                 self_type &m = (*this) ();
 752                 return m.find2 (1, index1 (), 0);
 753             }
 754             BOOST_UBLAS_INLINE
 755 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 756             typename self_type::
 757 #endif
 758             iterator2 end () const {
 759                 self_type &m = (*this) ();
 760                 return m.find2 (1, index1 (), m.size2 ());
 761             }
 762             BOOST_UBLAS_INLINE
 763 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 764             typename self_type::
 765 #endif
 766             reverse_iterator2 rbegin () const {
 767                 return reverse_iterator2 (end ());
 768             }
 769             BOOST_UBLAS_INLINE
 770 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 771             typename self_type::
 772 #endif
 773             reverse_iterator2 rend () const {
 774                 return reverse_iterator2 (begin ());
 775             }
 776 #endif
 777
 778             // Indices
 779             BOOST_UBLAS_INLINE
 780             size_type index1 () const {
 781                 BOOST_UBLAS_CHECK (*this != (*this) ().find1 (0, (*this) ().size1 (), j_), bad_index ());
 782                 if (rank_ == 1) {
 783                     BOOST_UBLAS_CHECK (layout_type::index_M (itv_.index (), it_.index ()) < (*this) ().size1 (), bad_index ());
 784                     return layout_type::index_M (itv_.index (), it_.index ());
 785                 } else {
 786                     return i_;
 787                 }
 788             }
 789             BOOST_UBLAS_INLINE
 790             size_type index2 () const {
 791                 BOOST_UBLAS_CHECK (*this != (*this) ().find1 (0, (*this) ().size1 (), j_), bad_index ());
 792                 if (rank_ == 1) {
 793                     BOOST_UBLAS_CHECK (layout_type::index_m (itv_.index (), it_.index ()) < (*this) ().size2 (), bad_index ());
 794                     return layout_type::index_m (itv_.index (), it_.index ());
 795                 } else {
 796                     return j_;
 797                 }
 798             }
 799
 800             // Assignment
 801             BOOST_UBLAS_INLINE
 802             iterator1 &operator = (const iterator1 &it) {
 803                 container_reference<self_type>::assign (&it ());
 804                 rank_ = it.rank_;
 805                 i_ = it.i_;
 806                 j_ = it.j_;
 807                 itv_ = it.itv_;
 808                 it_ = it.it_;
 809                 return *this;
 810             }
 811
 812             // Comparison
 813             BOOST_UBLAS_INLINE
 814             bool operator == (const iterator1 &it) const {
 815                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
 816                 // BOOST_UBLAS_CHECK (rank_ == it.rank_, internal_logic ());
 817                 if (rank_ == 1 || it.rank_ == 1) {
 818                     return it_ == it.it_;
 819                 } else {
 820                     return i_ == it.i_ && j_ == it.j_;
 821                 }
 822             }
 823
 824         private:
 825             int rank_;
 826             size_type i_;
 827             size_type j_;
 828             vectoriterator_type itv_;
 829             subiterator_type it_;
 830
 831             friend class const_iterator1;
 832         };
 833
 834         BOOST_UBLAS_INLINE
 835         iterator1 begin1 () {
 836             return find1 (0, 0, 0);
 837         }
 838         BOOST_UBLAS_INLINE
 839         iterator1 end1 () {
 840             return find1 (0, size1_, 0);
 841         }
 842
 843         class const_iterator2:
 844             public container_const_reference<generalized_vector_of_vector>,
 845             public bidirectional_iterator_base<sparse_bidirectional_iterator_tag,
 846                                                const_iterator2, value_type> {
 847         public:
 848             typedef typename generalized_vector_of_vector::difference_type difference_type;
 849             typedef typename generalized_vector_of_vector::value_type value_type;
 850             typedef typename generalized_vector_of_vector::const_reference reference;
 851             typedef const typename generalized_vector_of_vector::pointer pointer;
 852
 853             typedef const_iterator1 dual_iterator_type;
 854             typedef const_reverse_iterator1 dual_reverse_iterator_type;
 855
 856             // Construction and destruction
 857             BOOST_UBLAS_INLINE
 858             const_iterator2 ():
 859                 container_const_reference<self_type> (), rank_ (), i_ (), j_ (), itv_ (), it_ () {}
 860             BOOST_UBLAS_INLINE
 861             const_iterator2 (const self_type &m, int rank, size_type i, size_type j, const const_vectoriterator_type &itv, const const_subiterator_type &it):
 862                 container_const_reference<self_type> (m), rank_ (rank), i_ (i), j_ (j), itv_ (itv), it_ (it) {}
 863             BOOST_UBLAS_INLINE
 864             const_iterator2 (const iterator2 &it):
 865                 container_const_reference<self_type> (it ()), rank_ (it.rank_), i_ (it.i_), j_ (it.j_), itv_ (it.itv_), it_ (it.it_) {}
 866
 867             // Arithmetic
 868             BOOST_UBLAS_INLINE
 869             const_iterator2 &operator ++ () {
 870                 if (rank_ == 1 && layout_type::fast_j ())
 871                     ++ it_;
 872                 else {
 873                     const self_type &m = (*this) ();
 874                     j_ = index2 () + 1;
 875                     if (rank_ == 1 && ++ itv_ == m.end2 ().itv_)
 876                         *this = m.find2 (rank_, i_, j_, 1);
 877                     else if (rank_ == 1) {
 878                         it_ = (*itv_).begin ();
 879                         if (it_ == (*itv_).end () || index1 () != i_)
 880                             *this = m.find2 (rank_, i_, j_, 1);
 881                     }
 882                 }
 883                 return *this;
 884             }
 885             BOOST_UBLAS_INLINE
 886             const_iterator2 &operator -- () {
 887                 if (rank_ == 1 && layout_type::fast_j ())
 888                     -- it_;
 889                 else {
 890                     const self_type &m = (*this) ();
 891                     j_ = index2 () - 1;
 892                     if (rank_ == 1 && -- itv_ == m.end2 ().itv_)
 893                         *this = m.find2 (rank_, i_, j_, -1);
 894                     else if (rank_ == 1) {
 895                         it_ = (*itv_).begin ();
 896                         if (it_ == (*itv_).end () || index1 () != i_)
 897                             *this = m.find2 (rank_, i_, j_, -1);
 898                     }
 899                 }
 900                 return *this;
 901             }
 902
 903             // Dereference
 904             BOOST_UBLAS_INLINE
 905             const_reference operator * () const {
 906                 BOOST_UBLAS_CHECK (index1 () < (*this) ().size1 (), bad_index ());
 907                 BOOST_UBLAS_CHECK (index2 () < (*this) ().size2 (), bad_index ());
 908                 if (rank_ == 1) {
 909                     return *it_;
 910                 } else {
 911                     return (*this) () (i_, j_);
 912                 }
 913             }
 914
 915 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
 916             BOOST_UBLAS_INLINE
 917 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 918             typename self_type::
 919 #endif
 920             const_iterator1 begin () const {
 921                 const self_type &m = (*this) ();
 922                 return m.find1 (1, 0, index2 ());
 923             }
 924             BOOST_UBLAS_INLINE
 925 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 926             typename self_type::
 927 #endif
 928             const_iterator1 end () const {
 929                 const self_type &m = (*this) ();
 930                 return m.find1 (1, m.size1 (), index2 ());
 931             }
 932             BOOST_UBLAS_INLINE
 933 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 934             typename self_type::
 935 #endif
 936             const_reverse_iterator1 rbegin () const {
 937                 return const_reverse_iterator1 (end ());
 938             }
 939             BOOST_UBLAS_INLINE
 940 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
 941             typename self_type::
 942 #endif
 943             const_reverse_iterator1 rend () const {
 944                 return const_reverse_iterator1 (begin ());
 945             }
 946 #endif
 947
 948             // Indices
 949             BOOST_UBLAS_INLINE
 950             size_type index1 () const {
 951                 BOOST_UBLAS_CHECK (*this != (*this) ().find2 (0, i_, (*this) ().size2 ()), bad_index ());
 952                 if (rank_ == 1) {
 953                     BOOST_UBLAS_CHECK (layout_type::index_M (itv_.index (), it_.index ()) < (*this) ().size1 (), bad_index ());
 954                     return layout_type::index_M (itv_.index (), it_.index ());
 955                 } else {
 956                     return i_;
 957                 }
 958             }
 959             BOOST_UBLAS_INLINE
 960             size_type index2 () const {
 961                 BOOST_UBLAS_CHECK (*this != (*this) ().find2 (0, i_, (*this) ().size2 ()), bad_index ());
 962                 if (rank_ == 1) {
 963                     BOOST_UBLAS_CHECK (layout_type::index_m (itv_.index (), it_.index ()) < (*this) ().size2 (), bad_index ());
 964                     return layout_type::index_m (itv_.index (), it_.index ());
 965                 } else {
 966                     return j_;
 967                 }
 968             }
 969
 970             // Assignment
 971             BOOST_UBLAS_INLINE
 972             const_iterator2 &operator = (const const_iterator2 &it) {
 973                 container_const_reference<self_type>::assign (&it ());
 974                 rank_ = it.rank_;
 975                 i_ = it.i_;
 976                 j_ = it.j_;
 977                 itv_ = it.itv_;
 978                 it_ = it.it_;
 979                 return *this;
 980             }
 981
 982             // Comparison
 983             BOOST_UBLAS_INLINE
 984             bool operator == (const const_iterator2 &it) const {
 985                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
 986                 // BOOST_UBLAS_CHECK (rank_ == it.rank_, internal_logic ());
 987                 if (rank_ == 1 || it.rank_ == 1) {
 988                     return it_ == it.it_;
 989                 } else {
 990                     return i_ == it.i_ && j_ == it.j_;
 991                 }
 992             }
 993
 994         private:
 995             int rank_;
 996             size_type i_;
 997             size_type j_;
 998             const_vectoriterator_type itv_;
 999             const_subiterator_type it_;
1000         };
1001
1002         BOOST_UBLAS_INLINE
1003         const_iterator2 begin2 () const {
1004             return find2 (0, 0, 0);
1005         }
1006         BOOST_UBLAS_INLINE
1007         const_iterator2 end2 () const {
1008             return find2 (0, 0, size2_);
1009         }
1010
1011         class iterator2:
1012             public container_reference<generalized_vector_of_vector>,
1013             public bidirectional_iterator_base<sparse_bidirectional_iterator_tag,
1014                                                iterator2, value_type> {
1015         public:
1016             typedef typename generalized_vector_of_vector::difference_type difference_type;
1017             typedef typename generalized_vector_of_vector::value_type value_type;
1018             typedef typename generalized_vector_of_vector::true_reference reference;
1019             typedef typename generalized_vector_of_vector::pointer pointer;
1020
1021             typedef iterator1 dual_iterator_type;
1022             typedef reverse_iterator1 dual_reverse_iterator_type;
1023
1024             // Construction and destruction
1025             BOOST_UBLAS_INLINE
1026             iterator2 ():
1027                 container_reference<self_type> (), rank_ (), i_ (), j_ (), itv_ (), it_ () {}
1028             BOOST_UBLAS_INLINE
1029             iterator2 (self_type &m, int rank, size_type i, size_type j, const vectoriterator_type &itv, const subiterator_type &it):
1030                 container_reference<self_type> (m), rank_ (rank), i_ (i), j_ (j), itv_ (itv), it_ (it) {}
1031
1032             // Arithmetic
1033             BOOST_UBLAS_INLINE
1034             iterator2 &operator ++ () {
1035                 if (rank_ == 1 && layout_type::fast_j ())
1036                     ++ it_;
1037                 else {
1038                     self_type &m = (*this) ();
1039                     j_ = index2 () + 1;
1040                     if (rank_ == 1 && ++ itv_ == m.end2 ().itv_)
1041                         *this = m.find2 (rank_, i_, j_, 1);
1042                     else if (rank_ == 1) {
1043                         it_ = (*itv_).begin ();
1044                         if (it_ == (*itv_).end () || index1 () != i_)
1045                             *this = m.find2 (rank_, i_, j_, 1);
1046                     }
1047                 }
1048                 return *this;
1049             }
1050             BOOST_UBLAS_INLINE
1051             iterator2 &operator -- () {
1052                 if (rank_ == 1 && layout_type::fast_j ())
1053                     -- it_;
1054                 else {
1055                     self_type &m = (*this) ();
1056                     j_ = index2 () - 1;
1057                     if (rank_ == 1 && -- itv_ == m.end2 ().itv_)
1058                         *this = m.find2 (rank_, i_, j_, -1);
1059                     else if (rank_ == 1) {
1060                         it_ = (*itv_).begin ();
1061                         if (it_ == (*itv_).end () || index1 () != i_)
1062                             *this = m.find2 (rank_, i_, j_, -1);
1063                     }
1064                 }
1065                 return *this;
1066             }
1067
1068             // Dereference
1069             BOOST_UBLAS_INLINE
1070             true_reference operator * () const {
1071                 BOOST_UBLAS_CHECK (index1 () < (*this) ().size1 (), bad_index ());
1072                 BOOST_UBLAS_CHECK (index2 () < (*this) ().size2 (), bad_index ());
1073                 if (rank_ == 1) {
1074                     return *it_;
1075                 } else {
1076                     return (*this) ().at_element (i_, j_);
1077                 }
1078             }
1079
1080 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1081             BOOST_UBLAS_INLINE
1082 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1083             typename self_type::
1084 #endif
1085             iterator1 begin () const {
1086                 self_type &m = (*this) ();
1087                 return m.find1 (1, 0, index2 ());
1088             }
1089             BOOST_UBLAS_INLINE
1090 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1091             typename self_type::
1092 #endif
1093             iterator1 end () const {
1094                 self_type &m = (*this) ();
1095                 return m.find1 (1, m.size1 (), index2 ());
1096             }
1097             BOOST_UBLAS_INLINE
1098 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1099             typename self_type::
1100 #endif
1101             reverse_iterator1 rbegin () const {
1102                 return reverse_iterator1 (end ());
1103             }
1104             BOOST_UBLAS_INLINE
1105 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1106             typename self_type::
1107 #endif
1108             reverse_iterator1 rend () const {
1109                 return reverse_iterator1 (begin ());
1110             }
1111 #endif
1112
1113             // Indices
1114             BOOST_UBLAS_INLINE
1115             size_type index1 () const {
1116                 BOOST_UBLAS_CHECK (*this != (*this) ().find2 (0, i_, (*this) ().size2 ()), bad_index ());
1117                 if (rank_ == 1) {
1118                     BOOST_UBLAS_CHECK (layout_type::index_M (itv_.index (), it_.index ()) < (*this) ().size1 (), bad_index ());
1119                     return layout_type::index_M (itv_.index (), it_.index ());
1120                 } else {
1121                     return i_;
1122                 }
1123             }
1124             BOOST_UBLAS_INLINE
1125             size_type index2 () const {
1126                 BOOST_UBLAS_CHECK (*this != (*this) ().find2 (0, i_, (*this) ().size2 ()), bad_index ());
1127                 if (rank_ == 1) {
1128                     BOOST_UBLAS_CHECK (layout_type::index_m (itv_.index (), it_.index ()) < (*this) ().size2 (), bad_index ());
1129                     return layout_type::index_m (itv_.index (), it_.index ());
1130                 } else {
1131                     return j_;
1132                 }
1133             }
1134
1135             // Assignment
1136             BOOST_UBLAS_INLINE
1137             iterator2 &operator = (const iterator2 &it) {
1138                 container_reference<self_type>::assign (&it ());
1139                 rank_ = it.rank_;
1140                 i_ = it.i_;
1141                 j_ = it.j_;
1142                 itv_ = it.itv_;
1143                 it_ = it.it_;
1144                 return *this;
1145             }
1146
1147             // Comparison
1148             BOOST_UBLAS_INLINE
1149             bool operator == (const iterator2 &it) const {
1150                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1151                 // BOOST_UBLAS_CHECK (rank_ == it.rank_, internal_logic ());
1152                 if (rank_ == 1 || it.rank_ == 1) {
1153                     return it_ == it.it_;
1154                 } else {
1155                     return i_ == it.i_ && j_ == it.j_;
1156                 }
1157             }
1158
1159         private:
1160             int rank_;
1161             size_type i_;
1162             size_type j_;
1163             vectoriterator_type itv_;
1164             subiterator_type it_;
1165
1166             friend class const_iterator2;
1167         };
1168
1169         BOOST_UBLAS_INLINE
1170         iterator2 begin2 () {
1171             return find2 (0, 0, 0);
1172         }
1173         BOOST_UBLAS_INLINE
1174         iterator2 end2 () {
1175             return find2 (0, 0, size2_);
1176         }
1177
1178         // Reverse iterators
1179
1180         BOOST_UBLAS_INLINE
1181         const_reverse_iterator1 rbegin1 () const {
1182             return const_reverse_iterator1 (end1 ());
1183         }
1184         BOOST_UBLAS_INLINE
1185         const_reverse_iterator1 rend1 () const {
1186             return const_reverse_iterator1 (begin1 ());
1187         }
1188
1189         BOOST_UBLAS_INLINE
1190         reverse_iterator1 rbegin1 () {
1191             return reverse_iterator1 (end1 ());
1192         }
1193         BOOST_UBLAS_INLINE
1194         reverse_iterator1 rend1 () {
1195             return reverse_iterator1 (begin1 ());
1196         }
1197
1198         BOOST_UBLAS_INLINE
1199         const_reverse_iterator2 rbegin2 () const {
1200             return const_reverse_iterator2 (end2 ());
1201         }
1202         BOOST_UBLAS_INLINE
1203         const_reverse_iterator2 rend2 () const {
1204             return const_reverse_iterator2 (begin2 ());
1205         }
1206
1207         BOOST_UBLAS_INLINE
1208         reverse_iterator2 rbegin2 () {
1209             return reverse_iterator2 (end2 ());
1210         }
1211         BOOST_UBLAS_INLINE
1212         reverse_iterator2 rend2 () {
1213             return reverse_iterator2 (begin2 ());
1214         }
1215
1216          // Serialization
1217         template<class Archive>
1218         void serialize(Archive & ar, const unsigned int /* file_version */){
1219
1220             // we need to copy to a collection_size_type to get a portable
1221             // and efficient serialization
1222             serialization::collection_size_type s1 (size1_);
1223             serialization::collection_size_type s2 (size2_);
1224
1225             // serialize the sizes
1226             ar & serialization::make_nvp("size1",s1)
1227                & serialization::make_nvp("size2",s2);
1228
1229             // copy the values back if loading
1230             if (Archive::is_loading::value) {
1231                 size1_ = s1;
1232                 size2_ = s2;
1233             }
1234
1235             ar & serialization::make_nvp("data", data_);
1236
1237             storage_invariants();
1238         }
1239
1240     private:
1241         void storage_invariants () const
1242         {
1243             BOOST_UBLAS_CHECK (layout_type::size_M (size1_, size2_) + 1 == data_.size (), internal_logic ());
1244             BOOST_UBLAS_CHECK (data ().begin () != data ().end (), internal_logic ());
1245
1246         }
1247         size_type size1_;
1248         size_type size2_;
1249         array_type data_;
1250         static const value_type zero_;
1251     };
1252
1253     template<class T, class L, class A>
1254     const typename generalized_vector_of_vector<T, L, A>::value_type generalized_vector_of_vector<T, L, A>::zero_ = value_type/*zero*/();
1255
1256 }}}
1257
1258 #endif