boost/numeric/ublas/storage.hpp

   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_STORAGE_H
  14 #define BOOST_UBLAS_STORAGE_H
  15
  16 #include <algorithm>
  17 #ifdef BOOST_UBLAS_SHALLOW_ARRAY_ADAPTOR
  18 #include <boost/shared_array.hpp>
  19 #endif
  20
  21 #include <boost/serialization/array.hpp>
  22 #include <boost/serialization/collection_size_type.hpp>
  23 #include <boost/serialization/nvp.hpp>
  24
  25 #include <boost/numeric/ublas/exception.hpp>
  26 #include <boost/numeric/ublas/traits.hpp>
  27 #include <boost/numeric/ublas/detail/iterator.hpp>
  28
  29
  30 namespace boost { namespace numeric { namespace ublas {
  31
  32
  33     // Base class for Storage Arrays - see the Barton Nackman trick
  34     template<class E>
  35     class storage_array:
  36         private nonassignable {
  37     };
  38
  39
  40     // Unbounded array - with allocator
  41     template<class T, class ALLOC>
  42     class unbounded_array:
  43         public storage_array<unbounded_array<T, ALLOC> > {
  44
  45         typedef unbounded_array<T, ALLOC> self_type;
  46     public:
  47         typedef ALLOC allocator_type;
  48         typedef typename ALLOC::size_type size_type;
  49         typedef typename ALLOC::difference_type difference_type;
  50         typedef T value_type;
  51         typedef const T &const_reference;
  52         typedef T &reference;
  53         typedef const T *const_pointer;
  54         typedef T *pointer;
  55         typedef const_pointer const_iterator;
  56         typedef pointer iterator;
  57
  58         // Construction and destruction
  59         explicit BOOST_UBLAS_INLINE
  60         unbounded_array (const ALLOC &a = ALLOC()):
  61             alloc_ (a), size_ (0) {
  62             data_ = 0;
  63         }
  64         explicit BOOST_UBLAS_INLINE
  65         unbounded_array (size_type size, const ALLOC &a = ALLOC()):
  66             alloc_(a), size_ (size) {
  67           if (size_) {
  68               data_ = alloc_.allocate (size_);
  69               if (! detail::has_trivial_constructor<T>::value) {
  70                   for (pointer d = data_; d != data_ + size_; ++d)
  71                       alloc_.construct(d, value_type());
  72               }
  73           }
  74           else
  75               data_ = 0;
  76         }
  77         // No value initialised, but still be default constructed
  78         BOOST_UBLAS_INLINE
  79         unbounded_array (size_type size, const value_type &init, const ALLOC &a = ALLOC()):
  80             alloc_ (a), size_ (size) {
  81             if (size_) {
  82                 data_ = alloc_.allocate (size_);
  83                 std::uninitialized_fill (begin(), end(), init);
  84             }
  85             else
  86                 data_ = 0;
  87         }
  88         BOOST_UBLAS_INLINE
  89         unbounded_array (const unbounded_array &c):
  90             storage_array<unbounded_array<T, ALLOC> >(),
  91             alloc_ (c.alloc_), size_ (c.size_) {
  92             if (size_) {
  93                 data_ = alloc_.allocate (size_);
  94                 std::uninitialized_copy (c.begin(), c.end(), begin());
  95             }
  96             else
  97                 data_ = 0;
  98         }
  99         BOOST_UBLAS_INLINE
 100         ~unbounded_array () {
 101             if (size_) {
 102                 if (! detail::has_trivial_destructor<T>::value) {
 103                     // std::_Destroy (begin(), end(), alloc_);
 104                     const iterator i_end = end();
 105                     for (iterator i = begin (); i != i_end; ++i) {
 106                         iterator_destroy (i);
 107                     }
 108                 }
 109                 alloc_.deallocate (data_, size_);
 110             }
 111         }
 112
 113         // Resizing
 114     private:
 115         BOOST_UBLAS_INLINE
 116         void resize_internal (const size_type size, const value_type init, const bool preserve) {
 117             if (size != size_) {
 118                 pointer p_data = data_;
 119                 if (size) {
 120                     data_ = alloc_.allocate (size);
 121                     if (preserve) {
 122                         pointer si = p_data;
 123                         pointer di = data_;
 124                         if (size < size_) {
 125                             for (; di != data_ + size; ++di) {
 126                                 alloc_.construct (di, *si);
 127                                 ++si;
 128                             }
 129                         }
 130                         else {
 131                             for (pointer si = p_data; si != p_data + size_; ++si) {
 132                                 alloc_.construct (di, *si);
 133                                 ++di;
 134                             }
 135                             for (; di != data_ + size; ++di) {
 136                                 alloc_.construct (di, init);
 137                             }
 138                         }
 139                     }
 140                     else {
 141                         if (! detail::has_trivial_constructor<T>::value) {
 142                             for (pointer di = data_; di != data_ + size; ++di)
 143                                 alloc_.construct (di, value_type());
 144                         }
 145                     }
 146                 }
 147
 148                 if (size_) {
 149                     if (! detail::has_trivial_destructor<T>::value) {
 150                         for (pointer si = p_data; si != p_data + size_; ++si)
 151                             alloc_.destroy (si);
 152                     }
 153                     alloc_.deallocate (p_data, size_);
 154                 }
 155
 156                 if (!size)
 157                     data_ = 0;
 158                 size_ = size;
 159             }
 160         }
 161     public:
 162         BOOST_UBLAS_INLINE
 163         void resize (size_type size) {
 164             resize_internal (size, value_type (), false);
 165         }
 166         BOOST_UBLAS_INLINE
 167         void resize (size_type size, value_type init) {
 168             resize_internal (size, init, true);
 169         }
 170
 171         // Random Access Container
 172         BOOST_UBLAS_INLINE
 173         size_type max_size () const {
 174             return ALLOC ().max_size();
 175         }
 176
 177         BOOST_UBLAS_INLINE
 178         bool empty () const {
 179             return size_ == 0;
 180         }
 181
 182         BOOST_UBLAS_INLINE
 183         size_type size () const {
 184             return size_;
 185         }
 186
 187         // Element access
 188         BOOST_UBLAS_INLINE
 189         const_reference operator [] (size_type i) const {
 190             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 191             return data_ [i];
 192         }
 193         BOOST_UBLAS_INLINE
 194         reference operator [] (size_type i) {
 195             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 196             return data_ [i];
 197         }
 198
 199         // Assignment
 200         BOOST_UBLAS_INLINE
 201         unbounded_array &operator = (const unbounded_array &a) {
 202             if (this != &a) {
 203                 resize (a.size_);
 204                 std::copy (a.data_, a.data_ + a.size_, data_);
 205             }
 206             return *this;
 207         }
 208         BOOST_UBLAS_INLINE
 209         unbounded_array &assign_temporary (unbounded_array &a) {
 210             swap (a);
 211             return *this;
 212         }
 213
 214         // Swapping
 215         BOOST_UBLAS_INLINE
 216         void swap (unbounded_array &a) {
 217             if (this != &a) {
 218                 std::swap (size_, a.size_);
 219                 std::swap (data_, a.data_);
 220             }
 221         }
 222         BOOST_UBLAS_INLINE
 223         friend void swap (unbounded_array &a1, unbounded_array &a2) {
 224             a1.swap (a2);
 225         }
 226
 227         BOOST_UBLAS_INLINE
 228         const_iterator begin () const {
 229             return data_;
 230         }
 231         BOOST_UBLAS_INLINE
 232         const_iterator end () const {
 233             return data_ + size_;
 234         }
 235
 236         BOOST_UBLAS_INLINE
 237         iterator begin () {
 238             return data_;
 239         }
 240         BOOST_UBLAS_INLINE
 241         iterator end () {
 242             return data_ + size_;
 243         }
 244
 245         // Reverse iterators
 246         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
 247         typedef std::reverse_iterator<iterator> reverse_iterator;
 248
 249         BOOST_UBLAS_INLINE
 250         const_reverse_iterator rbegin () const {
 251             return const_reverse_iterator (end ());
 252         }
 253         BOOST_UBLAS_INLINE
 254         const_reverse_iterator rend () const {
 255             return const_reverse_iterator (begin ());
 256         }
 257         BOOST_UBLAS_INLINE
 258         reverse_iterator rbegin () {
 259             return reverse_iterator (end ());
 260         }
 261         BOOST_UBLAS_INLINE
 262         reverse_iterator rend () {
 263             return reverse_iterator (begin ());
 264         }
 265
 266         // Allocator
 267         allocator_type get_allocator () {
 268             return alloc_;
 269         }
 270
 271     private:
 272         friend class boost::serialization::access;
 273
 274         // Serialization
 275         template<class Archive>
 276         void serialize(Archive & ar, const unsigned int version)
 277         {
 278             serialization::collection_size_type s(size_);
 279             ar & serialization::make_nvp("size",s);
 280             if ( Archive::is_loading::value ) {
 281                 resize(s);
 282             }
 283             ar & serialization::make_array(data_, s);
 284         }
 285
 286     private:
 287         // Handle explict destroy on a (possibly indexed) iterator
 288         BOOST_UBLAS_INLINE
 289         static void iterator_destroy (iterator &i) {
 290             (&(*i)) -> ~value_type ();
 291         }
 292         ALLOC alloc_;
 293         size_type size_;
 294         pointer data_;
 295     };
 296
 297     // Bounded array - with allocator for size_type and difference_type
 298     template<class T, std::size_t N, class ALLOC>
 299     class bounded_array:
 300         public storage_array<bounded_array<T, N, ALLOC> > {
 301
 302         typedef bounded_array<T, N, ALLOC> self_type;
 303     public:
 304         // No allocator_type as ALLOC is not used for allocation
 305         typedef typename ALLOC::size_type size_type;
 306         typedef typename ALLOC::difference_type difference_type;
 307         typedef T value_type;
 308         typedef const T &const_reference;
 309         typedef T &reference;
 310         typedef const T *const_pointer;
 311         typedef T *pointer;
 312         typedef const_pointer const_iterator;
 313         typedef pointer iterator;
 314
 315         // Construction and destruction
 316         BOOST_UBLAS_INLINE
 317         bounded_array ():
 318             size_ (0) /*, data_ ()*/ {   // size 0 - use bounded_vector to default construct with size N
 319         }
 320         explicit BOOST_UBLAS_INLINE
 321         bounded_array (size_type size):
 322             size_ (size) /*, data_ ()*/ {
 323             BOOST_UBLAS_CHECK (size_ <= N, bad_size ());
 324             // data_ (an array) elements are already default constructed
 325         }
 326         BOOST_UBLAS_INLINE
 327         bounded_array (size_type size, const value_type &init):
 328             size_ (size) /*, data_ ()*/ {
 329             BOOST_UBLAS_CHECK (size_ <= N, bad_size ());
 330             // ISSUE elements should be value constructed here, but we must fill instead as already default constructed
 331             std::fill (begin(), end(), init) ;
 332         }
 333         BOOST_UBLAS_INLINE
 334         bounded_array (const bounded_array &c):
 335             size_ (c.size_)  {
 336             // ISSUE elements should be copy constructed here, but we must copy instead as already default constructed
 337             std::copy (c.begin(), c.end(), begin());
 338         }
 339
 340         // Resizing
 341         BOOST_UBLAS_INLINE
 342         void resize (size_type size) {
 343             BOOST_UBLAS_CHECK (size <= N, bad_size ());
 344             size_ = size;
 345         }
 346         BOOST_UBLAS_INLINE
 347         void resize (size_type size, value_type init) {
 348             BOOST_UBLAS_CHECK (size <= N, bad_size ());
 349             if (size > size_)
 350                 std::fill (data_ + size_, data_ + size, init);
 351             size_ = size;
 352         }
 353
 354         // Random Access Container
 355         BOOST_UBLAS_INLINE
 356         size_type max_size () const {
 357             return ALLOC ().max_size();
 358         }
 359
 360         BOOST_UBLAS_INLINE
 361         bool empty () const {
 362             return size_ == 0;
 363         }
 364
 365         BOOST_UBLAS_INLINE
 366         size_type size () const {
 367             return size_;
 368         }
 369
 370         // Element access
 371         BOOST_UBLAS_INLINE
 372         const_reference operator [] (size_type i) const {
 373             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 374             return data_ [i];
 375         }
 376         BOOST_UBLAS_INLINE
 377         reference operator [] (size_type i) {
 378             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 379             return data_ [i];
 380         }
 381
 382         // Assignment
 383         BOOST_UBLAS_INLINE
 384         bounded_array &operator = (const bounded_array &a) {
 385             if (this != &a) {
 386                 resize (a.size_);
 387                 std::copy (a.data_, a.data_ + a.size_, data_);
 388             }
 389             return *this;
 390         }
 391         BOOST_UBLAS_INLINE
 392         bounded_array &assign_temporary (bounded_array &a) {
 393             *this = a;
 394             return *this;
 395         }
 396
 397         // Swapping
 398         BOOST_UBLAS_INLINE
 399         void swap (bounded_array &a) {
 400             if (this != &a) {
 401                 std::swap (size_, a.size_);
 402                 std::swap_ranges (data_, data_ + (std::max) (size_, a.size_), a.data_);
 403             }
 404         }
 405         BOOST_UBLAS_INLINE
 406         friend void swap (bounded_array &a1, bounded_array &a2) {
 407             a1.swap (a2);
 408         }
 409
 410         BOOST_UBLAS_INLINE
 411         const_iterator begin () const {
 412             return data_;
 413         }
 414         BOOST_UBLAS_INLINE
 415         const_iterator end () const {
 416             return data_ + size_;
 417         }
 418
 419         BOOST_UBLAS_INLINE
 420         iterator begin () {
 421             return data_;
 422         }
 423         BOOST_UBLAS_INLINE
 424         iterator end () {
 425             return data_ + size_;
 426         }
 427
 428         // Reverse iterators
 429         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
 430         typedef std::reverse_iterator<iterator> reverse_iterator;
 431
 432         BOOST_UBLAS_INLINE
 433         const_reverse_iterator rbegin () const {
 434             return const_reverse_iterator (end ());
 435         }
 436         BOOST_UBLAS_INLINE
 437         const_reverse_iterator rend () const {
 438             return const_reverse_iterator (begin ());
 439         }
 440         BOOST_UBLAS_INLINE
 441         reverse_iterator rbegin () {
 442             return reverse_iterator (end ());
 443         }
 444         BOOST_UBLAS_INLINE
 445         reverse_iterator rend () {
 446             return reverse_iterator (begin ());
 447         }
 448
 449     private:
 450         // Serialization
 451         friend class boost::serialization::access;
 452
 453         template<class Archive>
 454         void serialize(Archive & ar, const unsigned int version)
 455         {
 456             serialization::collection_size_type s(size_);
 457             ar & serialization::make_nvp("size", s);
 458             if ( Archive::is_loading::value ) {
 459                 if (s > N) bad_size("too large size in bounded_array::load()\n").raise();
 460                 resize(s);
 461             }
 462             ar & serialization::make_array(data_, s);
 463         }
 464
 465     private:
 466         size_type size_;
 467         BOOST_UBLAS_BOUNDED_ARRAY_ALIGN value_type data_ [N];
 468     };
 469
 470
 471     // Array adaptor with normal deep copy semantics of elements
 472     template<class T>
 473     class array_adaptor:
 474         public storage_array<array_adaptor<T> > {
 475
 476         typedef array_adaptor<T> self_type;
 477     public:
 478         typedef std::size_t size_type;
 479         typedef std::ptrdiff_t difference_type;
 480         typedef T value_type;
 481         typedef const T &const_reference;
 482         typedef T &reference;
 483         typedef const T *const_pointer;
 484         typedef T *pointer;
 485
 486         // Construction and destruction
 487         BOOST_UBLAS_INLINE
 488         array_adaptor ():
 489             size_ (0), own_ (true), data_ (new value_type [0]) {
 490         }
 491         explicit BOOST_UBLAS_INLINE
 492         array_adaptor (size_type size):
 493             size_ (size), own_ (true), data_ (new value_type [size]) {
 494         }
 495         BOOST_UBLAS_INLINE
 496         array_adaptor (size_type size, const value_type &init):
 497             size_ (size), own_ (true), data_ (new value_type [size]) {
 498             std::fill (data_, data_ + size_, init);
 499         }
 500         BOOST_UBLAS_INLINE
 501         array_adaptor (size_type size, pointer data):
 502             size_ (size), own_ (false), data_ (data) {}
 503         BOOST_UBLAS_INLINE
 504         array_adaptor (const array_adaptor &a):
 505             storage_array<self_type> (),
 506             size_ (a.size_), own_ (true), data_ (new value_type [a.size_]) {
 507             *this = a;
 508         }
 509         BOOST_UBLAS_INLINE
 510         ~array_adaptor () {
 511             if (own_) {
 512                 delete [] data_;
 513             }
 514         }
 515
 516         // Resizing
 517     private:
 518         BOOST_UBLAS_INLINE
 519         void resize_internal (size_type size, value_type init, bool preserve = true) {
 520            if (size != size_) {
 521                 pointer data = new value_type [size];
 522                 if (preserve) {
 523                     std::copy (data_, data_ + (std::min) (size, size_), data);
 524                     std::fill (data + (std::min) (size, size_), data + size, init);
 525                 }
 526                 if (own_)
 527                     delete [] data_;
 528                 size_ = size;
 529                 own_ = true;
 530                 data_ = data;
 531             }
 532         }
 533         BOOST_UBLAS_INLINE
 534         void resize_internal (size_type size, pointer data, value_type init, bool preserve = true) {
 535             if (data != data_) {
 536                 if (preserve) {
 537                     std::copy (data_, data_ + (std::min) (size, size_), data);
 538                     std::fill (data + (std::min) (size, size_), data + size, init);
 539                 }
 540                 if (own_)
 541                     delete [] data_;
 542                 own_ = false;
 543                 data_ = data;
 544             }
 545             else {
 546                 std::fill (data + (std::min) (size, size_), data + size, init);
 547             }
 548             size_ = size;
 549         }
 550     public:
 551         BOOST_UBLAS_INLINE
 552         void resize (size_type size) {
 553             resize_internal (size, value_type (), false);
 554         }
 555         BOOST_UBLAS_INLINE
 556         void resize (size_type size, value_type init) {
 557             resize_internal (size, init, true);
 558         }
 559         BOOST_UBLAS_INLINE
 560         void resize (size_type size, pointer data) {
 561             resize_internal (size, data, value_type (), false);
 562         }
 563         BOOST_UBLAS_INLINE
 564         void resize (size_type size, pointer data, value_type init) {
 565             resize_internal (size, data, init, true);
 566         }
 567
 568         BOOST_UBLAS_INLINE
 569         size_type size () const {
 570             return size_;
 571         }
 572
 573         // Element access
 574         BOOST_UBLAS_INLINE
 575         const_reference operator [] (size_type i) const {
 576             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 577             return data_ [i];
 578         }
 579         BOOST_UBLAS_INLINE
 580         reference operator [] (size_type i) {
 581             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 582             return data_ [i];
 583         }
 584
 585         // Assignment
 586         BOOST_UBLAS_INLINE
 587         array_adaptor &operator = (const array_adaptor &a) {
 588             if (this != &a) {
 589                 resize (a.size_);
 590                 std::copy (a.data_, a.data_ + a.size_, data_);
 591             }
 592             return *this;
 593         }
 594         BOOST_UBLAS_INLINE
 595         array_adaptor &assign_temporary (array_adaptor &a) {
 596             if (own_ && a.own_)
 597                 swap (a);
 598             else
 599                 *this = a;
 600             return *this;
 601         }
 602
 603         // Swapping
 604         BOOST_UBLAS_INLINE
 605         void swap (array_adaptor &a) {
 606             if (this != &a) {
 607                 std::swap (size_, a.size_);
 608                 std::swap (own_, a.own_);
 609                 std::swap (data_, a.data_);
 610             }
 611         }
 612         BOOST_UBLAS_INLINE
 613         friend void swap (array_adaptor &a1, array_adaptor &a2) {
 614             a1.swap (a2);
 615         }
 616
 617         // Iterators simply are pointers.
 618
 619         typedef const_pointer const_iterator;
 620
 621         BOOST_UBLAS_INLINE
 622         const_iterator begin () const {
 623             return data_;
 624         }
 625         BOOST_UBLAS_INLINE
 626         const_iterator end () const {
 627             return data_ + size_;
 628         }
 629
 630         typedef pointer iterator;
 631
 632         BOOST_UBLAS_INLINE
 633         iterator begin () {
 634             return data_;
 635         }
 636         BOOST_UBLAS_INLINE
 637         iterator end () {
 638             return data_ + size_;
 639         }
 640
 641         // Reverse iterators
 642         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
 643         typedef std::reverse_iterator<iterator> reverse_iterator;
 644
 645         BOOST_UBLAS_INLINE
 646         const_reverse_iterator rbegin () const {
 647             return const_reverse_iterator (end ());
 648         }
 649         BOOST_UBLAS_INLINE
 650         const_reverse_iterator rend () const {
 651             return const_reverse_iterator (begin ());
 652         }
 653         BOOST_UBLAS_INLINE
 654         reverse_iterator rbegin () {
 655             return reverse_iterator (end ());
 656         }
 657         BOOST_UBLAS_INLINE
 658         reverse_iterator rend () {
 659             return reverse_iterator (begin ());
 660         }
 661
 662     private:
 663         size_type size_;
 664         bool own_;
 665         pointer data_;
 666     };
 667
 668 #ifdef BOOST_UBLAS_SHALLOW_ARRAY_ADAPTOR
 669     // Array adaptor with shallow (reference) copy semantics of elements.
 670     // shared_array is used to maintain reference counts.
 671     // This class breaks the normal copy semantics for a storage container and is very dangerous!
 672     template<class T>
 673     class shallow_array_adaptor:
 674         public storage_array<shallow_array_adaptor<T> > {
 675
 676         typedef shallow_array_adaptor<T> self_type;
 677
 678         template<class TT>
 679         struct leaker {
 680             typedef void result_type;
 681             typedef TT *argument_type;
 682
 683             BOOST_UBLAS_INLINE
 684             result_type operator () (argument_type x) {}
 685         };
 686
 687     public:
 688         typedef std::size_t size_type;
 689         typedef std::ptrdiff_t difference_type;
 690         typedef T value_type;
 691         typedef const T &const_reference;
 692         typedef T &reference;
 693         typedef const T *const_pointer;
 694         typedef T *pointer;
 695
 696         // Construction and destruction
 697         BOOST_UBLAS_INLINE
 698         shallow_array_adaptor ():
 699             size_ (0), own_ (true), data_ (new value_type [0]) {
 700         }
 701         explicit BOOST_UBLAS_INLINE
 702         shallow_array_adaptor (size_type size):
 703             size_ (size), own_ (true), data_ (new value_type [size]) {
 704         }
 705         BOOST_UBLAS_INLINE
 706         shallow_array_adaptor (size_type size, const value_type &init):
 707             size_ (size), own_ (true), data_ (new value_type [size]) {
 708             std::fill (data_.get (), data_.get () + size_, init);
 709         }
 710         BOOST_UBLAS_INLINE
 711         shallow_array_adaptor (size_type size, pointer data):
 712             size_ (size), own_ (false), data_ (data, leaker<value_type> ()) {}
 713
 714         BOOST_UBLAS_INLINE
 715         shallow_array_adaptor (const shallow_array_adaptor &a):
 716             storage_array<self_type> (),
 717             size_ (a.size_), own_ (a.own_), data_ (a.data_) {}
 718
 719         BOOST_UBLAS_INLINE
 720         ~shallow_array_adaptor () {
 721         }
 722
 723         // Resizing
 724     private:
 725         BOOST_UBLAS_INLINE
 726         void resize_internal (size_type size, value_type init, bool preserve = true) {
 727             if (size != size_) {
 728                 shared_array<value_type> data (new value_type [size]);
 729                 if (preserve) {
 730                     std::copy (data_.get (), data_.get () + (std::min) (size, size_), data.get ());
 731                     std::fill (data.get () + (std::min) (size, size_), data.get () + size, init);
 732                 }
 733                 size_ = size;
 734                 data_ = data;
 735             }
 736         }
 737         BOOST_UBLAS_INLINE
 738         void resize_internal (size_type size, pointer data, value_type init, bool preserve = true) {
 739             if (preserve) {
 740                 std::copy (data_.get (), data_.get () + (std::min) (size, size_), data);
 741                 std::fill (data + (std::min) (size, size_), data + size, init);
 742             }
 743             size_ = size;
 744             data_ = data;
 745         }
 746     public:
 747         BOOST_UBLAS_INLINE
 748         void resize (size_type size) {
 749             resize_internal (size, value_type (), false);
 750         }
 751         BOOST_UBLAS_INLINE
 752         void resize (size_type size, value_type init) {
 753             resize_internal (size, init, true);
 754         }
 755         BOOST_UBLAS_INLINE
 756         void resize (size_type size, pointer data) {
 757             resize_internal (size, data, value_type (), false);
 758         }
 759         BOOST_UBLAS_INLINE
 760         void resize (size_type size, pointer data, value_type init) {
 761             resize_internal (size, data, init, true);
 762         }
 763
 764         BOOST_UBLAS_INLINE
 765         size_type size () const {
 766             return size_;
 767         }
 768
 769         // Element access
 770         BOOST_UBLAS_INLINE
 771         const_reference operator [] (size_type i) const {
 772             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 773             return data_ [i];
 774         }
 775         BOOST_UBLAS_INLINE
 776         reference operator [] (size_type i) {
 777             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 778             return data_ [i];
 779         }
 780
 781         // Assignment
 782         BOOST_UBLAS_INLINE
 783         shallow_array_adaptor &operator = (const shallow_array_adaptor &a) {
 784             if (this != &a) {
 785                 resize (a.size_);
 786                 std::copy (a.data_.get (), a.data_.get () + a.size_, data_.get ());
 787             }
 788             return *this;
 789         }
 790         BOOST_UBLAS_INLINE
 791         shallow_array_adaptor &assign_temporary (shallow_array_adaptor &a) {
 792             if (own_ && a.own_)
 793                 swap (a);
 794             else
 795                 *this = a;
 796             return *this;
 797         }
 798
 799         // Swapping
 800         BOOST_UBLAS_INLINE
 801         void swap (shallow_array_adaptor &a) {
 802             if (this != &a) {
 803                 std::swap (size_, a.size_);
 804                 std::swap (own_, a.own_);
 805                 std::swap (data_, a.data_);
 806             }
 807         }
 808         BOOST_UBLAS_INLINE
 809         friend void swap (shallow_array_adaptor &a1, shallow_array_adaptor &a2) {
 810             a1.swap (a2);
 811         }
 812
 813         // Iterators simply are pointers.
 814
 815         typedef const_pointer const_iterator;
 816
 817         BOOST_UBLAS_INLINE
 818         const_iterator begin () const {
 819             return data_.get ();
 820         }
 821         BOOST_UBLAS_INLINE
 822         const_iterator end () const {
 823             return data_.get () + size_;
 824         }
 825
 826         typedef pointer iterator;
 827
 828         BOOST_UBLAS_INLINE
 829         iterator begin () {
 830             return data_.get ();
 831         }
 832         BOOST_UBLAS_INLINE
 833         iterator end () {
 834             return data_.get () + size_;
 835         }
 836
 837         // Reverse iterators
 838         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
 839         typedef std::reverse_iterator<iterator> reverse_iterator;
 840
 841         BOOST_UBLAS_INLINE
 842         const_reverse_iterator rbegin () const {
 843             return const_reverse_iterator (end ());
 844         }
 845         BOOST_UBLAS_INLINE
 846         const_reverse_iterator rend () const {
 847             return const_reverse_iterator (begin ());
 848         }
 849         BOOST_UBLAS_INLINE
 850         reverse_iterator rbegin () {
 851             return reverse_iterator (end ());
 852         }
 853         BOOST_UBLAS_INLINE
 854         reverse_iterator rend () {
 855             return reverse_iterator (begin ());
 856         }
 857
 858     private:
 859         size_type size_;
 860         bool own_;
 861         shared_array<value_type> data_;
 862     };
 863
 864 #endif
 865
 866
 867     // Range class
 868     template <class Z, class D>
 869     class basic_range {
 870         typedef basic_range<Z, D> self_type;
 871     public:
 872         typedef Z size_type;
 873         typedef D difference_type;
 874         typedef size_type value_type;
 875         typedef value_type const_reference;
 876         typedef const_reference reference;
 877         typedef const value_type *const_pointer;
 878         typedef value_type *pointer;
 879
 880         // Construction and destruction
 881         BOOST_UBLAS_INLINE
 882         basic_range ():
 883             start_ (0), size_ (0) {}
 884         BOOST_UBLAS_INLINE
 885         basic_range (size_type start, size_type stop):
 886             start_ (start), size_ (stop - start) {
 887             BOOST_UBLAS_CHECK (start_ <= stop, bad_index ());
 888         }
 889
 890         BOOST_UBLAS_INLINE
 891         size_type start () const {
 892             return start_;
 893         }
 894         BOOST_UBLAS_INLINE
 895         size_type size () const {
 896             return size_;
 897         }
 898
 899         // Random Access Container
 900         BOOST_UBLAS_INLINE
 901         size_type max_size () const {
 902             return size_;
 903         }
 904
 905         BOOST_UBLAS_INLINE
 906         bool empty () const {
 907             return size_ == 0;
 908         }
 909
 910         // Element access
 911         BOOST_UBLAS_INLINE
 912         const_reference operator () (size_type i) const {
 913             BOOST_UBLAS_CHECK (i < size_, bad_index ());
 914             return start_ + i;
 915         }
 916
 917         // Composition
 918         BOOST_UBLAS_INLINE
 919         basic_range compose (const basic_range &r) const {
 920             return basic_range (start_ + r.start_, start_ + r.start_ + r.size_);
 921         }
 922
 923         // Comparison
 924         BOOST_UBLAS_INLINE
 925         bool operator == (const basic_range &r) const {
 926             return start_ == r.start_ && size_ == r.size_;
 927         }
 928         BOOST_UBLAS_INLINE
 929         bool operator != (const basic_range &r) const {
 930             return ! (*this == r);
 931         }
 932
 933         // Iterator types
 934     private:
 935         // Use and index
 936         typedef size_type const_subiterator_type;
 937
 938     public:
 939 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
 940         typedef indexed_const_iterator<self_type, std::random_access_iterator_tag> const_iterator;
 941 #else
 942         class const_iterator:
 943             public container_const_reference<basic_range>,
 944             public random_access_iterator_base<std::random_access_iterator_tag,
 945                                                const_iterator, value_type> {
 946         public:
 947             typedef typename basic_range::value_type value_type;
 948             typedef typename basic_range::difference_type difference_type;
 949             typedef typename basic_range::const_reference reference;
 950             typedef typename basic_range::const_pointer pointer;
 951
 952             // Construction and destruction
 953             BOOST_UBLAS_INLINE
 954             const_iterator ():
 955                 container_const_reference<basic_range> (), it_ () {}
 956             BOOST_UBLAS_INLINE
 957             const_iterator (const basic_range &r, const const_subiterator_type &it):
 958                 container_const_reference<basic_range> (r), it_ (it) {}
 959
 960             // Arithmetic
 961             BOOST_UBLAS_INLINE
 962             const_iterator &operator ++ () {
 963                 ++ it_;
 964                 return *this;
 965             }
 966             BOOST_UBLAS_INLINE
 967             const_iterator &operator -- () {
 968                 BOOST_UBLAS_CHECK (it_ > 0, bad_index ());
 969                 -- it_;
 970                 return *this;
 971             }
 972             BOOST_UBLAS_INLINE
 973             const_iterator &operator += (difference_type n) {
 974                 BOOST_UBLAS_CHECK (n >= 0 || it_ >= size_type(-n), bad_index ());
 975                 it_ += n;
 976                 return *this;
 977             }
 978             BOOST_UBLAS_INLINE
 979             const_iterator &operator -= (difference_type n) {
 980                 BOOST_UBLAS_CHECK (n <= 0 || it_ >= size_type(n), bad_index ());
 981                 it_ -= n;
 982                 return *this;
 983             }
 984             BOOST_UBLAS_INLINE
 985             difference_type operator - (const const_iterator &it) const {
 986                 return it_ - it.it_;
 987             }
 988
 989             // Dereference
 990             BOOST_UBLAS_INLINE
 991             const_reference operator * () const {
 992                 BOOST_UBLAS_CHECK ((*this) ().start () <= it_, bad_index ());
 993                 BOOST_UBLAS_CHECK (it_ < (*this) ().start () + (*this) ().size (), bad_index ());
 994                 return it_;
 995             }
 996
 997             BOOST_UBLAS_INLINE
 998             const_reference operator [] (difference_type n) const {
 999                 return *(*this + n);
1000             }
1001
1002             // Index
1003             BOOST_UBLAS_INLINE
1004             size_type index () const {
1005                 BOOST_UBLAS_CHECK ((*this) ().start () <= it_, bad_index ());
1006                 BOOST_UBLAS_CHECK (it_ < (*this) ().start () + (*this) ().size (), bad_index ());
1007                 return it_ - (*this) ().start ();
1008             }
1009
1010             // Assignment
1011             BOOST_UBLAS_INLINE
1012             const_iterator &operator = (const const_iterator &it) {
1013                 // Comeau recommends...
1014                 this->assign (&it ());
1015                 it_ = it.it_;
1016                 return *this;
1017             }
1018
1019             // Comparison
1020             BOOST_UBLAS_INLINE
1021             bool operator == (const const_iterator &it) const {
1022                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1023                 return it_ == it.it_;
1024             }
1025             BOOST_UBLAS_INLINE
1026             bool operator < (const const_iterator &it) const {
1027                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1028                 return it_ < it.it_;
1029             }
1030
1031         private:
1032             const_subiterator_type it_;
1033         };
1034 #endif
1035
1036         BOOST_UBLAS_INLINE
1037         const_iterator begin () const {
1038             return const_iterator (*this, start_);
1039         }
1040         BOOST_UBLAS_INLINE
1041         const_iterator end () const {
1042             return const_iterator (*this, start_ + size_);
1043         }
1044
1045         // Reverse iterator
1046         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
1047
1048         BOOST_UBLAS_INLINE
1049         const_reverse_iterator rbegin () const {
1050             return const_reverse_iterator (end ());
1051         }
1052         BOOST_UBLAS_INLINE
1053         const_reverse_iterator rend () const {
1054             return const_reverse_iterator (begin ());
1055         }
1056
1057         BOOST_UBLAS_INLINE
1058         basic_range preprocess (size_type size) const {
1059             if (this != &all_)
1060                 return *this;
1061             return basic_range (0, size);
1062         }
1063         static
1064         BOOST_UBLAS_INLINE
1065         const basic_range &all () {
1066             return all_;
1067         }
1068
1069     private:
1070         size_type start_;
1071         size_type size_;
1072         static const basic_range all_;
1073     };
1074
1075     template <class Z, class D>
1076     const basic_range<Z,D> basic_range<Z,D>::all_  (0, size_type (-1));
1077
1078
1079     // Slice class
1080     template <class Z, class D>
1081     class basic_slice {
1082         typedef basic_slice<Z, D> self_type;
1083     public:
1084         typedef Z size_type;
1085         typedef D difference_type;
1086         typedef size_type value_type;
1087         typedef value_type const_reference;
1088         typedef const_reference reference;
1089         typedef const value_type *const_pointer;
1090         typedef value_type *pointer;
1091
1092         // Construction and destruction
1093         BOOST_UBLAS_INLINE
1094         basic_slice ():
1095             start_ (0), stride_ (0), size_ (0) {}
1096         BOOST_UBLAS_INLINE
1097         basic_slice (size_type start, difference_type stride, size_type size):
1098             start_ (start), stride_ (stride), size_ (size) {}
1099
1100         BOOST_UBLAS_INLINE
1101         size_type start () const {
1102             return start_;
1103         }
1104         BOOST_UBLAS_INLINE
1105         difference_type stride () const {
1106             return stride_;
1107         }
1108         BOOST_UBLAS_INLINE
1109         size_type size () const {
1110             return size_;
1111         }
1112
1113         // Random Access Container
1114         BOOST_UBLAS_INLINE
1115         size_type max_size () const {
1116             return size_;
1117         }
1118
1119         BOOST_UBLAS_INLINE
1120         bool empty () const {
1121             return size_ == 0;
1122         }
1123
1124         // Element access
1125         BOOST_UBLAS_INLINE
1126         const_reference operator () (size_type i) const {
1127             BOOST_UBLAS_CHECK (i < size_, bad_index ());
1128             BOOST_UBLAS_CHECK (stride_ >= 0 || start_ >= i * -stride_, bad_index ());
1129             return start_ + i * stride_;
1130         }
1131
1132         // Composition
1133         BOOST_UBLAS_INLINE
1134         basic_slice compose (const basic_range<size_type, difference_type> &r) const {
1135             BOOST_UBLAS_CHECK (stride_ >=0 || start_ >= -stride_ * r.start(), bad_index ());
1136             return basic_slice (start_ + stride_ * r.start (), stride_, r.size ());
1137         }
1138         BOOST_UBLAS_INLINE
1139         basic_slice compose (const basic_slice &s) const {
1140             BOOST_UBLAS_CHECK (stride_ >=0 || start_ >= -stride_ * s.start_, bad_index ());
1141             return basic_slice (start_ + stride_ * s.start_, stride_ * s.stride_, s.size_);
1142         }
1143
1144         // Comparison
1145         BOOST_UBLAS_INLINE
1146         bool operator == (const basic_slice &s) const {
1147             return start_ == s.start_ && stride_ == s.stride_ && size_ == s.size_;
1148         }
1149         BOOST_UBLAS_INLINE
1150         bool operator != (const basic_slice &s) const {
1151             return ! (*this == s);
1152         }
1153
1154         // Iterator types
1155     private:
1156         // Use and index
1157         typedef size_type const_subiterator_type;
1158
1159     public:
1160 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1161         typedef indexed_const_iterator<self_type, std::random_access_iterator_tag> const_iterator;
1162 #else
1163         class const_iterator:
1164             public container_const_reference<basic_slice>,
1165             public random_access_iterator_base<std::random_access_iterator_tag,
1166                                                const_iterator, value_type> {
1167         public:
1168             typedef typename basic_slice::value_type value_type;
1169             typedef typename basic_slice::difference_type difference_type;
1170             typedef typename basic_slice::const_reference reference;
1171             typedef typename basic_slice::const_pointer pointer;
1172
1173             // Construction and destruction
1174             BOOST_UBLAS_INLINE
1175             const_iterator ():
1176                 container_const_reference<basic_slice> (), it_ () {}
1177             BOOST_UBLAS_INLINE
1178             const_iterator (const basic_slice &s, const const_subiterator_type &it):
1179                 container_const_reference<basic_slice> (s), it_ (it) {}
1180
1181             // Arithmetic
1182             BOOST_UBLAS_INLINE
1183             const_iterator &operator ++ () {
1184                 ++it_;
1185                 return *this;
1186             }
1187             BOOST_UBLAS_INLINE
1188             const_iterator &operator -- () {
1189                 BOOST_UBLAS_CHECK (it_ > 0, bad_index ());
1190                 --it_;
1191                 return *this;
1192             }
1193             BOOST_UBLAS_INLINE
1194             const_iterator &operator += (difference_type n) {
1195                 BOOST_UBLAS_CHECK (n >= 0 || it_ >= size_type(-n), bad_index ());
1196                 it_ += n;
1197                 return *this;
1198             }
1199             BOOST_UBLAS_INLINE
1200             const_iterator &operator -= (difference_type n) {
1201                 BOOST_UBLAS_CHECK (n <= 0 || it_ >= size_type(n), bad_index ());
1202                 it_ -= n;
1203                 return *this;
1204             }
1205             BOOST_UBLAS_INLINE
1206             difference_type operator - (const const_iterator &it) const {
1207                 return it_ - it.it_;
1208             }
1209
1210             // Dereference
1211             BOOST_UBLAS_INLINE
1212             const_reference operator * () const {
1213                 BOOST_UBLAS_CHECK (it_ < (*this) ().size (), bad_index ());
1214                 return (*this) ().start () + it_* (*this) ().stride ();
1215             }
1216
1217             BOOST_UBLAS_INLINE
1218             const_reference operator [] (difference_type n) const {
1219                 return *(*this + n);
1220             }
1221
1222             // Index
1223             BOOST_UBLAS_INLINE
1224             size_type index () const {
1225                 BOOST_UBLAS_CHECK (it_ < (*this) ().size (), bad_index ());
1226                 return it_;
1227             }
1228
1229             // Assignment
1230             BOOST_UBLAS_INLINE
1231             const_iterator &operator = (const const_iterator &it) {
1232                 // Comeau recommends...
1233                 this->assign (&it ());
1234                 it_ = it.it_;
1235                 return *this;
1236             }
1237
1238             // Comparison
1239             BOOST_UBLAS_INLINE
1240             bool operator == (const const_iterator &it) const {
1241                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1242                 return it_ == it.it_;
1243             }
1244             BOOST_UBLAS_INLINE
1245             bool operator < (const const_iterator &it) const {
1246                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1247                 return it_ < it.it_;
1248             }
1249
1250         private:
1251             const_subiterator_type it_;
1252         };
1253 #endif
1254
1255         BOOST_UBLAS_INLINE
1256         const_iterator begin () const {
1257             return const_iterator (*this, 0);
1258         }
1259         BOOST_UBLAS_INLINE
1260         const_iterator end () const {
1261             return const_iterator (*this, size_);
1262         }
1263
1264         // Reverse iterator
1265         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
1266
1267         BOOST_UBLAS_INLINE
1268         const_reverse_iterator rbegin () const {
1269             return const_reverse_iterator (end ());
1270         }
1271         BOOST_UBLAS_INLINE
1272         const_reverse_iterator rend () const {
1273             return const_reverse_iterator (begin ());
1274         }
1275
1276         BOOST_UBLAS_INLINE
1277         basic_slice preprocess (size_type size) const {
1278             if (this != &all_)
1279                 return *this;
1280             return basic_slice (0, 1, size);
1281         }
1282         static
1283         BOOST_UBLAS_INLINE
1284         const basic_slice &all () {
1285             return all_;
1286         }
1287
1288     private:
1289         size_type start_;
1290         difference_type stride_;
1291         size_type size_;
1292         static const basic_slice all_;
1293     };
1294
1295     template <class Z, class D>
1296     const basic_slice<Z,D> basic_slice<Z,D>::all_  (0, 1, size_type (-1));
1297
1298
1299     // Indirect array class
1300     template<class A>
1301     class indirect_array {
1302         typedef indirect_array<A> self_type;
1303     public:
1304         typedef A array_type;
1305         typedef const A const_array_type;
1306         typedef typename A::size_type size_type;
1307         typedef typename A::difference_type difference_type;
1308         typedef typename A::value_type value_type;
1309         typedef typename A::const_reference const_reference;
1310         typedef typename A::reference reference;
1311         typedef typename A::const_pointer const_pointer;
1312         typedef typename A::pointer pointer;
1313
1314         // Construction and destruction
1315         BOOST_UBLAS_INLINE
1316         indirect_array ():
1317             size_ (), data_ () {}
1318         explicit BOOST_UBLAS_INLINE
1319         indirect_array (size_type size):
1320             size_ (size), data_ (size) {}
1321         BOOST_UBLAS_INLINE
1322         indirect_array (size_type size, const array_type &data):
1323             size_ (size), data_ (data) {}
1324         BOOST_UBLAS_INLINE
1325         indirect_array (pointer start, pointer stop):
1326             size_ (stop - start), data_ (stop - start) {
1327             std::copy (start, stop, data_.begin ());
1328         }
1329
1330         BOOST_UBLAS_INLINE
1331         size_type size () const {
1332             return size_;
1333         }
1334         BOOST_UBLAS_INLINE
1335         const_array_type data () const {
1336             return data_;
1337         }
1338         BOOST_UBLAS_INLINE
1339         array_type data () {
1340             return data_;
1341         }
1342
1343         // Random Access Container
1344         BOOST_UBLAS_INLINE
1345         size_type max_size () const {
1346             return size_;
1347         }
1348
1349         BOOST_UBLAS_INLINE
1350         bool empty () const {
1351             return data_.size () == 0;
1352         }
1353
1354         // Element access
1355         BOOST_UBLAS_INLINE
1356         const_reference operator () (size_type i) const {
1357             BOOST_UBLAS_CHECK (i < size_, bad_index ());
1358             return data_ [i];
1359         }
1360         BOOST_UBLAS_INLINE
1361         reference operator () (size_type i) {
1362             BOOST_UBLAS_CHECK (i < size_, bad_index ());
1363             return data_ [i];
1364         }
1365
1366         BOOST_UBLAS_INLINE
1367         const_reference operator [] (size_type i) const {
1368             return (*this) (i);
1369         }
1370         BOOST_UBLAS_INLINE
1371         reference operator [] (size_type i) {
1372             return (*this) (i);
1373         }
1374
1375         // Composition
1376         BOOST_UBLAS_INLINE
1377         indirect_array compose (const basic_range<size_type, difference_type> &r) const {
1378             BOOST_UBLAS_CHECK (r.start () + r.size () <= size_, bad_size ());
1379             array_type data (r.size ());
1380             for (size_type i = 0; i < r.size (); ++ i)
1381                 data [i] = data_ [r.start () + i];
1382             return indirect_array (r.size (), data);
1383         }
1384         BOOST_UBLAS_INLINE
1385         indirect_array compose (const basic_slice<size_type, difference_type> &s) const {
1386             BOOST_UBLAS_CHECK (s.start () + s.stride () * (s.size () - (s.size () > 0)) <= size (), bad_size ());
1387             array_type data (s.size ());
1388             for (size_type i = 0; i < s.size (); ++ i)
1389                 data [i] = data_ [s.start () + s.stride () * i];
1390             return indirect_array (s.size (), data);
1391         }
1392         BOOST_UBLAS_INLINE
1393         indirect_array compose (const indirect_array &ia) const {
1394             array_type data (ia.size_);
1395             for (size_type i = 0; i < ia.size_; ++ i) {
1396                 BOOST_UBLAS_CHECK (ia.data_ [i] <= size_, bad_size ());
1397                 data [i] = data_ [ia.data_ [i]];
1398             }
1399             return indirect_array (ia.size_, data);
1400         }
1401
1402         // Comparison
1403         template<class OA>
1404         BOOST_UBLAS_INLINE
1405         bool operator == (const indirect_array<OA> &ia) const {
1406             if (size_ != ia.size_)
1407                 return false;
1408             for (size_type i = 0; i < BOOST_UBLAS_SAME (size_, ia.size_); ++ i)
1409                 if (data_ [i] != ia.data_ [i])
1410                     return false;
1411             return true;
1412         }
1413         template<class OA>
1414         BOOST_UBLAS_INLINE
1415         bool operator != (const indirect_array<OA> &ia) const {
1416             return ! (*this == ia);
1417         }
1418
1419         // Iterator types
1420     private:
1421         // Use a index difference
1422         typedef difference_type const_subiterator_type;
1423
1424     public:
1425 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1426         typedef indexed_const_iterator<indirect_array, std::random_access_iterator_tag> const_iterator;
1427 #else
1428         class const_iterator:
1429             public container_const_reference<indirect_array>,
1430             public random_access_iterator_base<std::random_access_iterator_tag,
1431                                                const_iterator, value_type> {
1432         public:
1433             typedef typename indirect_array::value_type value_type;
1434             typedef typename indirect_array::difference_type difference_type;
1435             typedef typename indirect_array::const_reference reference;
1436             typedef typename indirect_array::const_pointer pointer;
1437
1438             // Construction and destruction
1439             BOOST_UBLAS_INLINE
1440             const_iterator ():
1441                 container_const_reference<indirect_array> (), it_ () {}
1442             BOOST_UBLAS_INLINE
1443             const_iterator (const indirect_array &ia, const const_subiterator_type &it):
1444                 container_const_reference<indirect_array> (ia), it_ (it) {}
1445
1446             // Arithmetic
1447             BOOST_UBLAS_INLINE
1448             const_iterator &operator ++ () {
1449                 ++ it_;
1450                 return *this;
1451             }
1452             BOOST_UBLAS_INLINE
1453             const_iterator &operator -- () {
1454                 -- it_;
1455                 return *this;
1456             }
1457             BOOST_UBLAS_INLINE
1458             const_iterator &operator += (difference_type n) {
1459                 it_ += n;
1460                 return *this;
1461             }
1462             BOOST_UBLAS_INLINE
1463             const_iterator &operator -= (difference_type n) {
1464                 it_ -= n;
1465                 return *this;
1466             }
1467             BOOST_UBLAS_INLINE
1468             difference_type operator - (const const_iterator &it) const {
1469                 return it_ - it.it_;
1470             }
1471
1472             // Dereference
1473             BOOST_UBLAS_INLINE
1474             const_reference operator * () const {
1475                 return (*this) () (it_);
1476             }
1477
1478             BOOST_UBLAS_INLINE
1479             const_reference operator [] (difference_type n) const {
1480                 return *(*this + n);
1481             }
1482
1483             // Index
1484             BOOST_UBLAS_INLINE
1485             size_type index () const {
1486                 return it_;
1487             }
1488
1489             // Assignment
1490             BOOST_UBLAS_INLINE
1491             const_iterator &operator = (const const_iterator &it) {
1492                 // Comeau recommends...
1493                 this->assign (&it ());
1494                 it_ = it.it_;
1495                 return *this;
1496             }
1497
1498             // Comparison
1499             BOOST_UBLAS_INLINE
1500             bool operator == (const const_iterator &it) const {
1501                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1502                 return it_ == it.it_;
1503             }
1504             BOOST_UBLAS_INLINE
1505             bool operator < (const const_iterator &it) const {
1506                 BOOST_UBLAS_CHECK ((*this) () == it (), external_logic ());
1507                 return it_ < it.it_;
1508             }
1509
1510         private:
1511             const_subiterator_type it_;
1512         };
1513 #endif
1514
1515         BOOST_UBLAS_INLINE
1516         const_iterator begin () const {
1517             return const_iterator (*this, 0);
1518         }
1519         BOOST_UBLAS_INLINE
1520         const_iterator end () const {
1521             return const_iterator (*this, size_);
1522         }
1523
1524         // Reverse iterator
1525         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
1526
1527         BOOST_UBLAS_INLINE
1528         const_reverse_iterator rbegin () const {
1529             return const_reverse_iterator (end ());
1530         }
1531         BOOST_UBLAS_INLINE
1532         const_reverse_iterator rend () const {
1533             return const_reverse_iterator (begin ());
1534         }
1535
1536         BOOST_UBLAS_INLINE
1537         indirect_array preprocess (size_type size) const {
1538             if (this != &all_)
1539                 return *this;
1540             indirect_array ia (size);
1541             for (size_type i = 0; i < size; ++ i)
1542                ia (i) = i;
1543             return ia;
1544         }
1545         static
1546         BOOST_UBLAS_INLINE
1547         const indirect_array &all () {
1548             return all_;
1549         }
1550
1551     private:
1552         size_type size_;
1553         array_type data_;
1554         static const indirect_array all_;
1555     };
1556
1557     template<class A>
1558     const indirect_array<A> indirect_array<A>::all_;
1559
1560
1561
1562     // Gunter Winkler contributed the classes index_pair, index_pair_array,
1563     // index_triple and index_triple_array to enable inplace sort of parallel arrays.
1564
1565     template <class V>
1566     class index_pair :
1567         private boost::noncopyable,
1568         public container_reference<V> {
1569
1570         typedef index_pair<V> self_type;
1571     public:
1572         typedef typename V::size_type size_type;
1573
1574         BOOST_UBLAS_INLINE
1575         index_pair(V& v, size_type i) :
1576             container_reference<V>(v), i_(i),
1577             v1_(v.data1_[i]), v2_(v.data2_[i]),
1578             dirty_(false), is_copy_(false) {}
1579          BOOST_UBLAS_INLINE
1580         index_pair(const self_type& rhs) :
1581             container_reference<V>(rhs()), i_(0),
1582             v1_(rhs.v1_), v2_(rhs.v2_),
1583             dirty_(false), is_copy_(true) {}
1584          BOOST_UBLAS_INLINE
1585         ~index_pair() {
1586             if (dirty_ && (!is_copy_) ) {
1587                 (*this)().data1_[i_] = v1_;
1588                 (*this)().data2_[i_] = v2_;
1589             }
1590         }
1591
1592         BOOST_UBLAS_INLINE
1593         self_type& operator=(const self_type& rhs) {
1594             v1_ = rhs.v1_;
1595             v2_ = rhs.v2_;
1596             dirty_ = true;
1597             return *this;
1598         }
1599
1600         BOOST_UBLAS_INLINE
1601         void swap(self_type rhs) {
1602             self_type tmp(rhs);
1603             rhs = *this;
1604             *this = tmp;
1605         }
1606         BOOST_UBLAS_INLINE
1607         friend void swap(self_type lhs, self_type rhs) {
1608             lhs.swap(rhs);
1609         }
1610
1611         BOOST_UBLAS_INLINE
1612         bool equal(const self_type& rhs) const {
1613             return (v1_ == rhs.v1_);
1614         }
1615         BOOST_UBLAS_INLINE
1616         bool less(const self_type& rhs) const {
1617             return (v1_ < rhs.v1_);
1618         }
1619         BOOST_UBLAS_INLINE
1620         friend bool operator == (const self_type& lhs, const self_type& rhs) {
1621             return lhs.equal(rhs);
1622         }
1623         BOOST_UBLAS_INLINE
1624         friend bool operator != (const self_type& lhs, const self_type& rhs) {
1625             return !lhs.equal(rhs);
1626         }
1627         BOOST_UBLAS_INLINE
1628         friend bool operator < (const self_type& lhs, const self_type& rhs) {
1629             return lhs.less(rhs);
1630         }
1631         BOOST_UBLAS_INLINE
1632         friend bool operator >= (const self_type& lhs, const self_type& rhs) {
1633             return !lhs.less(rhs);
1634         }
1635         BOOST_UBLAS_INLINE
1636         friend bool operator > (const self_type& lhs, const self_type& rhs) {
1637             return rhs.less(lhs);
1638         }
1639         BOOST_UBLAS_INLINE
1640         friend bool operator <= (const self_type& lhs, const self_type& rhs) {
1641             return !rhs.less(lhs);
1642         }
1643
1644     private:
1645         size_type i_;
1646         typename V::value1_type v1_;
1647         typename V::value2_type v2_;
1648         bool dirty_;
1649         bool is_copy_;
1650      };
1651
1652     template <class V1, class V2>
1653     class index_pair_array:
1654         private boost::noncopyable {
1655
1656         typedef index_pair_array<V1, V2> self_type;
1657     public:
1658         typedef typename V1::value_type value1_type;
1659         typedef typename V2::value_type value2_type;
1660
1661         typedef typename V1::size_type size_type;
1662         typedef typename V1::difference_type difference_type;
1663         typedef index_pair<self_type> value_type;
1664         // There is nothing that can be referenced directly. Always return a copy of the index_pair
1665         typedef value_type reference;
1666         typedef const value_type const_reference;
1667
1668         BOOST_UBLAS_INLINE
1669         index_pair_array(size_type size, V1& data1, V2& data2) :
1670               size_(size),data1_(data1),data2_(data2) {}
1671
1672         BOOST_UBLAS_INLINE
1673         size_type size() const {
1674             return size_;
1675         }
1676
1677         BOOST_UBLAS_INLINE
1678         const_reference operator () (size_type i) const {
1679             return value_type((*this), i);
1680         }
1681         BOOST_UBLAS_INLINE
1682         reference operator () (size_type i) {
1683             return value_type((*this), i);
1684         }
1685
1686         typedef indexed_iterator<self_type, std::random_access_iterator_tag> iterator;
1687         typedef indexed_const_iterator<self_type, std::random_access_iterator_tag> const_iterator;
1688
1689         BOOST_UBLAS_INLINE
1690         iterator begin() {
1691             return iterator( (*this), 0);
1692         }
1693         BOOST_UBLAS_INLINE
1694         iterator end() {
1695             return iterator( (*this), size());
1696         }
1697
1698         BOOST_UBLAS_INLINE
1699         const_iterator begin() const {
1700             return const_iterator( (*this), 0);
1701         }
1702         BOOST_UBLAS_INLINE
1703         const_iterator end() const {
1704             return const_iterator( (*this), size());
1705         }
1706
1707         // unnecessary function:
1708         BOOST_UBLAS_INLINE
1709         bool equal(size_type i1, size_type i2) const {
1710             return data1_[i1] == data1_[i2];
1711         }
1712         BOOST_UBLAS_INLINE
1713         bool less(size_type i1, size_type i2) const {
1714             return data1_[i1] < data1_[i2];
1715         }
1716
1717         // gives a large speedup
1718         BOOST_UBLAS_INLINE
1719         friend void iter_swap(const iterator& lhs, const iterator& rhs) {
1720             const size_type i1 = lhs.index();
1721             const size_type i2 = rhs.index();
1722             std::swap(lhs().data1_[i1], rhs().data1_[i2]);
1723             std::swap(lhs().data2_[i1], rhs().data2_[i2]);
1724         }
1725
1726     private:
1727         size_type size_;
1728         V1& data1_;
1729         V2& data2_;
1730
1731         // friend class value_type;
1732         friend class index_pair<self_type>;
1733     };
1734
1735     template <class M>
1736     class index_triple :
1737         private boost::noncopyable,
1738         public container_reference<M> {
1739
1740         typedef index_triple<M> self_type;
1741     public:
1742         typedef typename M::size_type size_type;
1743
1744         BOOST_UBLAS_INLINE
1745         index_triple(M& m, size_type i) :
1746             container_reference<M>(m), i_(i),
1747             v1_(m.data1_[i]), v2_(m.data2_[i]), v3_(m.data3_[i]),
1748             dirty_(false), is_copy_(false) {}
1749         BOOST_UBLAS_INLINE
1750         index_triple(const self_type& rhs) :
1751             container_reference<M>(rhs()), i_(0),
1752             v1_(rhs.v1_), v2_(rhs.v2_), v3_(rhs.v3_),
1753             dirty_(false), is_copy_(true) {}
1754         BOOST_UBLAS_INLINE
1755         ~index_triple() {
1756             if (dirty_ && (!is_copy_) ) {
1757                 (*this)().data1_[i_] = v1_;
1758                 (*this)().data2_[i_] = v2_;
1759                 (*this)().data3_[i_] = v3_;
1760             }
1761         }
1762
1763         BOOST_UBLAS_INLINE
1764         self_type& operator=(const self_type& rhs) {
1765             v1_ = rhs.v1_;
1766             v2_ = rhs.v2_;
1767             v3_ = rhs.v3_;
1768             dirty_ = true;
1769             return *this;
1770         }
1771
1772         BOOST_UBLAS_INLINE
1773         void swap(self_type rhs) {
1774             self_type tmp(rhs);
1775             rhs = *this;
1776             *this = tmp;
1777         }
1778         BOOST_UBLAS_INLINE
1779         friend void swap(self_type lhs, self_type rhs) {
1780             lhs.swap(rhs);
1781         }
1782
1783         BOOST_UBLAS_INLINE
1784         bool equal(const self_type& rhs) const {
1785             return ((v1_ == rhs.v1_) && (v2_ == rhs.v2_));
1786         }
1787         BOOST_UBLAS_INLINE
1788         bool less(const self_type& rhs) const {
1789             return ((v1_ < rhs.v1_) ||
1790                     (v1_ == rhs.v1_ && v2_ < rhs.v2_));
1791         }
1792         BOOST_UBLAS_INLINE
1793         friend bool operator == (const self_type& lhs, const self_type& rhs) {
1794             return lhs.equal(rhs);
1795         }
1796         BOOST_UBLAS_INLINE
1797         friend bool operator != (const self_type& lhs, const self_type& rhs) {
1798             return !lhs.equal(rhs);
1799         }
1800         BOOST_UBLAS_INLINE
1801         friend bool operator < (const self_type& lhs, const self_type& rhs) {
1802             return lhs.less(rhs);
1803         }
1804         BOOST_UBLAS_INLINE
1805         friend bool operator >= (const self_type& lhs, const self_type& rhs) {
1806             return !lhs.less(rhs);
1807         }
1808         BOOST_UBLAS_INLINE
1809         friend bool operator > (const self_type& lhs, const self_type& rhs) {
1810             return rhs.less(lhs);
1811         }
1812         BOOST_UBLAS_INLINE
1813         friend bool operator <= (const self_type& lhs, const self_type& rhs) {
1814             return !rhs.less(lhs);
1815         }
1816
1817     private:
1818         size_type i_;
1819         typename M::value1_type v1_;
1820         typename M::value2_type v2_;
1821         typename M::value3_type v3_;
1822         bool dirty_;
1823         bool is_copy_;
1824     };
1825
1826     template <class V1, class V2, class V3>
1827     class index_triple_array:
1828         private boost::noncopyable {
1829
1830         typedef index_triple_array<V1, V2, V3> self_type;
1831     public:
1832         typedef typename V1::value_type value1_type;
1833         typedef typename V2::value_type value2_type;
1834         typedef typename V3::value_type value3_type;
1835
1836         typedef typename V1::size_type size_type;
1837         typedef typename V1::difference_type difference_type;
1838         typedef index_triple<self_type> value_type;
1839         // There is nothing that can be referenced directly. Always return a copy of the index_triple
1840         typedef value_type reference;
1841         typedef const value_type const_reference;
1842
1843         BOOST_UBLAS_INLINE
1844         index_triple_array(size_type size, V1& data1, V2& data2, V3& data3) :
1845               size_(size),data1_(data1),data2_(data2),data3_(data3) {}
1846
1847         BOOST_UBLAS_INLINE
1848         size_type size() const {
1849             return size_;
1850         }
1851
1852         BOOST_UBLAS_INLINE
1853         const_reference operator () (size_type i) const {
1854             return value_type((*this), i);
1855         }
1856         BOOST_UBLAS_INLINE
1857         reference operator () (size_type i) {
1858             return value_type((*this), i);
1859         }
1860
1861         typedef indexed_iterator<self_type, std::random_access_iterator_tag> iterator;
1862         typedef indexed_const_iterator<self_type, std::random_access_iterator_tag> const_iterator;
1863
1864         BOOST_UBLAS_INLINE
1865         iterator begin() {
1866             return iterator( (*this), 0);
1867         }
1868         BOOST_UBLAS_INLINE
1869         iterator end() {
1870             return iterator( (*this), size());
1871         }
1872
1873         BOOST_UBLAS_INLINE
1874         const_iterator begin() const {
1875             return const_iterator( (*this), 0);
1876         }
1877         BOOST_UBLAS_INLINE
1878         const_iterator end() const {
1879             return const_iterator( (*this), size());
1880         }
1881
1882         // unnecessary function:
1883         BOOST_UBLAS_INLINE
1884         bool equal(size_type i1, size_type i2) const {
1885             return ((data1_[i1] == data1_[i2]) && (data2_[i1] == data2_[i2]));
1886         }
1887         BOOST_UBLAS_INLINE
1888         bool less(size_type i1, size_type i2) const {
1889             return ((data1_[i1] < data1_[i2]) ||
1890                     (data1_[i1] == data1_[i2] && data2_[i1] < data2_[i2]));
1891         }
1892
1893         // gives a large speedup
1894         BOOST_UBLAS_INLINE
1895         friend void iter_swap(const iterator& lhs, const iterator& rhs) {
1896             const size_type i1 = lhs.index();
1897             const size_type i2 = rhs.index();
1898             std::swap(lhs().data1_[i1], rhs().data1_[i2]);
1899             std::swap(lhs().data2_[i1], rhs().data2_[i2]);
1900             std::swap(lhs().data3_[i1], rhs().data3_[i2]);
1901         }
1902
1903     private:
1904         size_type size_;
1905         V1& data1_;
1906         V2& data2_;
1907         V3& data3_;
1908
1909         // friend class value_type;
1910         friend class index_triple<self_type>;
1911     };
1912
1913 }}}
1914
1915 #endif