1 //
2 //  Copyright (c) 2000-2002
3 //  Joerg Walter, Mathias Koch
4 //
5 //  Permission to use, copy, modify, distribute and sell this software
6 //  and its documentation for any purpose is hereby granted without fee,
7 //  provided that the above copyright notice appear in all copies and
8 //  that both that copyright notice and this permission notice appear
9 //  in supporting documentation.  The authors make no representations
10 //  about the suitability of this software for any purpose.
11 //  It is provided "as is" without express or implied warranty.
12 //
13 //  The authors gratefully acknowledge the support of
14 //  GeNeSys mbH & Co. KG in producing this work.
15 //
16 
17 #ifndef _BOOST_UBLAS_MATRIX_PROXY_
18 #define _BOOST_UBLAS_MATRIX_PROXY_
19 
20 #include <boost/numeric/ublas/matrix_expression.hpp>
21 #include <boost/numeric/ublas/detail/vector_assign.hpp>
22 #include <boost/numeric/ublas/detail/matrix_assign.hpp>
23 #include <boost/numeric/ublas/detail/temporary.hpp>
24 
25 // Iterators based on ideas of Jeremy Siek
26 
27 namespace boost { namespace numeric { namespace ublas {
28 
29     // Matrix based row vector class
30     template<class M>
31     class matrix_row:
32         public vector_expression<matrix_row<M> > {
33 
34         typedef matrix_row<M> self_type;
35     public:
36 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
37         using vector_expression<self_type>::operator ();
38 #endif
39         typedef M matrix_type;
40         typedef typename M::size_type size_type;
41         typedef typename M::difference_type difference_type;
42         typedef typename M::value_type value_type;
43         typedef typename M::const_reference const_reference;
44         typedef typename boost::mpl::if_<boost::is_const<M>,
45                                           typename M::const_reference,
46                                           typename M::reference>::type reference;
47         typedef typename boost::mpl::if_<boost::is_const<M>,
48                                           typename M::const_closure_type,
49                                           typename M::closure_type>::type matrix_closure_type;
50         typedef const self_type const_closure_type;
51         typedef self_type closure_type;
52         typedef typename storage_restrict_traits<typename M::storage_category,
53                                                  dense_proxy_tag>::storage_category storage_category;
54 
55         // Construction and destruction
56         BOOST_UBLAS_INLINE
matrix_row(matrix_type & data,size_type i)57         matrix_row (matrix_type &data, size_type i):
58             data_ (data), i_ (i) {
59             // Early checking of preconditions here.
60             // BOOST_UBLAS_CHECK (i_ < data_.size1 (), bad_index ());
61         }
62 
63         // Accessors
64         BOOST_UBLAS_INLINE
size() const65         size_type size () const {
66             return data_.size2 ();
67         }
68         BOOST_UBLAS_INLINE
index() const69         size_type index () const {
70             return i_;
71         }
72 
73         // Storage accessors
74         BOOST_UBLAS_INLINE
data() const75         const matrix_closure_type &data () const {
76             return data_;
77         }
78         BOOST_UBLAS_INLINE
data()79         matrix_closure_type &data () {
80             return data_;
81         }
82 
83         // Element access
84 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
85         BOOST_UBLAS_INLINE
operator ()(size_type j) const86         const_reference operator () (size_type j) const {
87             return data_ (i_, j);
88         }
89         BOOST_UBLAS_INLINE
operator ()(size_type j)90         reference operator () (size_type j) {
91             return data_ (i_, j);
92         }
93 
94         BOOST_UBLAS_INLINE
operator [](size_type j) const95         const_reference operator [] (size_type j) const {
96             return (*this) (j);
97         }
98         BOOST_UBLAS_INLINE
operator [](size_type j)99         reference operator [] (size_type j) {
100             return (*this) (j);
101         }
102 #else
103         BOOST_UBLAS_INLINE
operator ()(size_type j) const104         reference operator () (size_type j) const {
105             return data_ (i_, j);
106         }
107 
108         BOOST_UBLAS_INLINE
operator [](size_type j) const109         reference operator [] (size_type j) const {
110             return (*this) (j);
111         }
112 #endif
113 
114         // Assignment
115         BOOST_UBLAS_INLINE
operator =(const matrix_row & mr)116         matrix_row &operator = (const matrix_row &mr) {
117             // ISSUE need a temporary, proxy can be overlaping alias
118             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (mr));
119             return *this;
120         }
121         BOOST_UBLAS_INLINE
assign_temporary(matrix_row & mr)122         matrix_row &assign_temporary (matrix_row &mr) {
123             // assign elements, proxied container remains the same
124             vector_assign<scalar_assign> (*this, mr);
125             return *this;
126         }
127         template<class AE>
128         BOOST_UBLAS_INLINE
operator =(const vector_expression<AE> & ae)129         matrix_row &operator = (const vector_expression<AE> &ae) {
130             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (ae));
131             return *this;
132         }
133         template<class AE>
134         BOOST_UBLAS_INLINE
assign(const vector_expression<AE> & ae)135         matrix_row &assign (const vector_expression<AE> &ae) {
136             vector_assign<scalar_assign> (*this, ae);
137             return *this;
138         }
139         template<class AE>
140         BOOST_UBLAS_INLINE
operator +=(const vector_expression<AE> & ae)141         matrix_row &operator += (const vector_expression<AE> &ae) {
142             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this + ae));
143             return *this;
144         }
145         template<class AE>
146         BOOST_UBLAS_INLINE
plus_assign(const vector_expression<AE> & ae)147         matrix_row &plus_assign (const vector_expression<AE> &ae) {
148             vector_assign<scalar_plus_assign> (*this, ae);
149             return *this;
150         }
151         template<class AE>
152         BOOST_UBLAS_INLINE
operator -=(const vector_expression<AE> & ae)153         matrix_row &operator -= (const vector_expression<AE> &ae) {
154             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this - ae));
155             return *this;
156         }
157         template<class AE>
158         BOOST_UBLAS_INLINE
minus_assign(const vector_expression<AE> & ae)159         matrix_row &minus_assign (const vector_expression<AE> &ae) {
160             vector_assign<scalar_minus_assign> (*this, ae);
161             return *this;
162         }
163         template<class AT>
164         BOOST_UBLAS_INLINE
operator *=(const AT & at)165         matrix_row &operator *= (const AT &at) {
166             vector_assign_scalar<scalar_multiplies_assign> (*this, at);
167             return *this;
168         }
169         template<class AT>
170         BOOST_UBLAS_INLINE
operator /=(const AT & at)171         matrix_row &operator /= (const AT &at) {
172             vector_assign_scalar<scalar_divides_assign> (*this, at);
173             return *this;
174         }
175 
176         // Closure comparison
177         BOOST_UBLAS_INLINE
same_closure(const matrix_row & mr) const178         bool same_closure (const matrix_row &mr) const {
179             return (*this).data_.same_closure (mr.data_);
180         }
181 
182         // Comparison
183         BOOST_UBLAS_INLINE
operator ==(const matrix_row & mr) const184         bool operator == (const matrix_row &mr) const {
185             return (*this).data_ == mr.data_ && index () == mr.index ();
186         }
187 
188         // Swapping
189         BOOST_UBLAS_INLINE
swap(matrix_row mr)190         void swap (matrix_row mr) {
191             if (this != &mr) {
192                 BOOST_UBLAS_CHECK (size () == mr.size (), bad_size ());
193                 // Sparse ranges may be nonconformant now.
194                 // std::swap_ranges (begin (), end (), mr.begin ());
195                 vector_swap<scalar_swap> (*this, mr);
196             }
197         }
198         BOOST_UBLAS_INLINE
swap(matrix_row mr1,matrix_row mr2)199         friend void swap (matrix_row mr1, matrix_row mr2) {
200             mr1.swap (mr2);
201         }
202 
203         // Iterator types
204     private:
205         typedef typename M::const_iterator2 const_subiterator_type;
206         typedef typename boost::mpl::if_<boost::is_const<M>,
207                                           typename M::const_iterator2,
208                                           typename M::iterator2>::type subiterator_type;
209 
210     public:
211 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
212         typedef indexed_iterator<matrix_row<matrix_type>,
213                                  typename subiterator_type::iterator_category> iterator;
214         typedef indexed_const_iterator<matrix_row<matrix_type>,
215                                        typename const_subiterator_type::iterator_category> const_iterator;
216 #else
217         class const_iterator;
218         class iterator;
219 #endif
220 
221         // Element lookup
222         BOOST_UBLAS_INLINE
find(size_type j) const223         const_iterator find (size_type j) const {
224             const_subiterator_type it2 (data_.find2 (1, i_, j));
225 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
226             return const_iterator (*this, it2.index2 ());
227 #else
228             return const_iterator (*this, it2);
229 #endif
230         }
231         BOOST_UBLAS_INLINE
find(size_type j)232         iterator find (size_type j) {
233             subiterator_type it2 (data_.find2 (1, i_, j));
234 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
235             return iterator (*this, it2.index2 ());
236 #else
237             return iterator (*this, it2);
238 #endif
239         }
240 
241 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
242         class const_iterator:
243             public container_const_reference<matrix_row>,
244             public iterator_base_traits<typename const_subiterator_type::iterator_category>::template
245                         iterator_base<const_iterator, value_type>::type {
246         public:
247             typedef typename const_subiterator_type::value_type value_type;
248             typedef typename const_subiterator_type::difference_type difference_type;
249             typedef typename const_subiterator_type::reference reference;
250             typedef typename const_subiterator_type::pointer pointer;
251 
252             // Construction and destruction
253             BOOST_UBLAS_INLINE
const_iterator()254             const_iterator ():
255                 container_const_reference<self_type> (), it_ () {}
256             BOOST_UBLAS_INLINE
const_iterator(const self_type & mr,const const_subiterator_type & it)257             const_iterator (const self_type &mr, const const_subiterator_type &it):
258                 container_const_reference<self_type> (mr), it_ (it) {}
259             BOOST_UBLAS_INLINE
const_iterator(const iterator & it)260             const_iterator (const iterator &it):
261                 container_const_reference<self_type> (it ()), it_ (it.it_) {}
262 
263             // Arithmetic
264             BOOST_UBLAS_INLINE
operator ++()265             const_iterator &operator ++ () {
266                 ++ it_;
267                 return *this;
268             }
269             BOOST_UBLAS_INLINE
operator --()270             const_iterator &operator -- () {
271                 -- it_;
272                 return *this;
273             }
274             BOOST_UBLAS_INLINE
operator +=(difference_type n)275             const_iterator &operator += (difference_type n) {
276                 it_ += n;
277                 return *this;
278             }
279             BOOST_UBLAS_INLINE
operator -=(difference_type n)280             const_iterator &operator -= (difference_type n) {
281                 it_ -= n;
282                 return *this;
283             }
284             BOOST_UBLAS_INLINE
operator -(const const_iterator & it) const285             difference_type operator - (const const_iterator &it) const {
286                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
287                 return it_ - it.it_;
288             }
289 
290             // Dereference
291             BOOST_UBLAS_INLINE
operator *() const292             const_reference operator * () const {
293                 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ());
294                 return *it_;
295             }
296 
297             // Index
298             BOOST_UBLAS_INLINE
index() const299             size_type index () const {
300                 return it_.index2 ();
301             }
302 
303             // Assignment
304             BOOST_UBLAS_INLINE
operator =(const const_iterator & it)305             const_iterator &operator = (const const_iterator &it) {
306                 container_const_reference<self_type>::assign (&it ());
307                 it_ = it.it_;
308                 return *this;
309             }
310 
311             // Comparison
312             BOOST_UBLAS_INLINE
operator ==(const const_iterator & it) const313             bool operator == (const const_iterator &it) const {
314                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
315                 return it_ == it.it_;
316             }
317             BOOST_UBLAS_INLINE
operator <(const const_iterator & it) const318             bool operator < (const const_iterator &it) const {
319                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
320                 return it_ < it.it_;
321             }
322 
323         private:
324             const_subiterator_type it_;
325         };
326 #endif
327 
328         BOOST_UBLAS_INLINE
begin() const329         const_iterator begin () const {
330             return find (0);
331         }
332         BOOST_UBLAS_INLINE
end() const333         const_iterator end () const {
334             return find (size ());
335         }
336 
337 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
338         class iterator:
339             public container_reference<matrix_row>,
340             public iterator_base_traits<typename subiterator_type::iterator_category>::template
341                         iterator_base<iterator, value_type>::type {
342         public:
343             typedef typename subiterator_type::value_type value_type;
344             typedef typename subiterator_type::difference_type difference_type;
345             typedef typename subiterator_type::reference reference;
346             typedef typename subiterator_type::pointer pointer;
347 
348             // Construction and destruction
349             BOOST_UBLAS_INLINE
iterator()350             iterator ():
351                 container_reference<self_type> (), it_ () {}
352             BOOST_UBLAS_INLINE
iterator(self_type & mr,const subiterator_type & it)353             iterator (self_type &mr, const subiterator_type &it):
354                 container_reference<self_type> (mr), it_ (it) {}
355 
356             // Arithmetic
357             BOOST_UBLAS_INLINE
operator ++()358             iterator &operator ++ () {
359                 ++ it_;
360                 return *this;
361             }
362             BOOST_UBLAS_INLINE
operator --()363             iterator &operator -- () {
364                 -- it_;
365                 return *this;
366             }
367             BOOST_UBLAS_INLINE
operator +=(difference_type n)368             iterator &operator += (difference_type n) {
369                 it_ += n;
370                 return *this;
371             }
372             BOOST_UBLAS_INLINE
operator -=(difference_type n)373             iterator &operator -= (difference_type n) {
374                 it_ -= n;
375                 return *this;
376             }
377             BOOST_UBLAS_INLINE
operator -(const iterator & it) const378             difference_type operator - (const iterator &it) const {
379                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
380                 return it_ - it.it_;
381             }
382 
383             // Dereference
384             BOOST_UBLAS_INLINE
operator *() const385             reference operator * () const {
386                 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ());
387                 return *it_;
388             }
389 
390             // Index
391             BOOST_UBLAS_INLINE
index() const392             size_type index () const {
393                 return it_.index2 ();
394             }
395 
396             // Assignment
397             BOOST_UBLAS_INLINE
operator =(const iterator & it)398             iterator &operator = (const iterator &it) {
399                 container_reference<self_type>::assign (&it ());
400                 it_ = it.it_;
401                 return *this;
402             }
403 
404             // Comparison
405             BOOST_UBLAS_INLINE
operator ==(const iterator & it) const406             bool operator == (const iterator &it) const {
407                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
408                 return it_ == it.it_;
409             }
410             BOOST_UBLAS_INLINE
operator <(const iterator & it) const411             bool operator < (const iterator &it) const {
412                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
413                 return it_ < it.it_;
414             }
415 
416         private:
417             subiterator_type it_;
418 
419             friend class const_iterator;
420         };
421 #endif
422 
423         BOOST_UBLAS_INLINE
begin()424         iterator begin () {
425             return find (0);
426         }
427         BOOST_UBLAS_INLINE
end()428         iterator end () {
429             return find (size ());
430         }
431 
432         // Reverse iterator
433         typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
434         typedef reverse_iterator_base<iterator> reverse_iterator;
435 
436         BOOST_UBLAS_INLINE
rbegin() const437         const_reverse_iterator rbegin () const {
438             return const_reverse_iterator (end ());
439         }
440         BOOST_UBLAS_INLINE
rend() const441         const_reverse_iterator rend () const {
442             return const_reverse_iterator (begin ());
443         }
444         BOOST_UBLAS_INLINE
rbegin()445         reverse_iterator rbegin () {
446             return reverse_iterator (end ());
447         }
448         BOOST_UBLAS_INLINE
rend()449         reverse_iterator rend () {
450             return reverse_iterator (begin ());
451         }
452 
453     private:
454         matrix_closure_type data_;
455         size_type i_;
456     };
457 
458     // Projections
459     template<class M>
460     BOOST_UBLAS_INLINE
row(M & data,typename M::size_type i)461     matrix_row<M> row (M &data, typename M::size_type i) {
462         return matrix_row<M> (data, i);
463     }
464     template<class M>
465     BOOST_UBLAS_INLINE
row(const M & data,typename M::size_type i)466     const matrix_row<const M> row (const M &data, typename M::size_type i) {
467         return matrix_row<const M> (data, i);
468     }
469 
470     // Specialize temporary
471     template <class M>
472     struct vector_temporary_traits< matrix_row<M> >
473     : vector_temporary_traits< M > {} ;
474 
475     // Matrix based column vector class
476     template<class M>
477     class matrix_column:
478         public vector_expression<matrix_column<M> > {
479 
480         typedef matrix_column<M> self_type;
481     public:
482 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
483         using vector_expression<self_type>::operator ();
484 #endif
485         typedef M matrix_type;
486         typedef typename M::size_type size_type;
487         typedef typename M::difference_type difference_type;
488         typedef typename M::value_type value_type;
489         typedef typename M::const_reference const_reference;
490         typedef typename boost::mpl::if_<boost::is_const<M>,
491                                           typename M::const_reference,
492                                           typename M::reference>::type reference;
493         typedef typename boost::mpl::if_<boost::is_const<M>,
494                                           typename M::const_closure_type,
495                                           typename M::closure_type>::type matrix_closure_type;
496         typedef const self_type const_closure_type;
497         typedef self_type closure_type;
498         typedef typename storage_restrict_traits<typename M::storage_category,
499                                                  dense_proxy_tag>::storage_category storage_category;
500 
501         // Construction and destruction
502         BOOST_UBLAS_INLINE
matrix_column(matrix_type & data,size_type j)503         matrix_column (matrix_type &data, size_type j):
504             data_ (data), j_ (j) {
505             // Early checking of preconditions here.
506             // BOOST_UBLAS_CHECK (j_ < data_.size2 (), bad_index ());
507         }
508 
509         // Accessors
510         BOOST_UBLAS_INLINE
size() const511         size_type size () const {
512             return data_.size1 ();
513         }
514         BOOST_UBLAS_INLINE
index() const515         size_type index () const {
516             return j_;
517         }
518 
519         // Storage accessors
520         BOOST_UBLAS_INLINE
data() const521         const matrix_closure_type &data () const {
522             return data_;
523         }
524         BOOST_UBLAS_INLINE
data()525         matrix_closure_type &data () {
526             return data_;
527         }
528 
529         // Element access
530 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
531         BOOST_UBLAS_INLINE
operator ()(size_type i) const532         const_reference operator () (size_type i) const {
533             return data_ (i, j_);
534         }
535         BOOST_UBLAS_INLINE
operator ()(size_type i)536         reference operator () (size_type i) {
537             return data_ (i, j_);
538         }
539 
540         BOOST_UBLAS_INLINE
operator [](size_type i) const541         const_reference operator [] (size_type i) const {
542             return (*this) (i);
543         }
544         BOOST_UBLAS_INLINE
operator [](size_type i)545         reference operator [] (size_type i) {
546             return (*this) (i);
547         }
548 #else
549         BOOST_UBLAS_INLINE
operator ()(size_type i) const550         reference operator () (size_type i) const {
551             return data_ (i, j_);
552         }
553 
554         BOOST_UBLAS_INLINE
operator [](size_type i) const555         reference operator [] (size_type i) const {
556             return (*this) (i);
557         }
558 #endif
559 
560         // Assignment
561         BOOST_UBLAS_INLINE
operator =(const matrix_column & mc)562         matrix_column &operator = (const matrix_column &mc) {
563             // ISSUE need a temporary, proxy can be overlaping alias
564             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (mc));
565             return *this;
566         }
567         BOOST_UBLAS_INLINE
assign_temporary(matrix_column & mc)568         matrix_column &assign_temporary (matrix_column &mc) {
569             // assign elements, proxied container remains the same
570             vector_assign<scalar_assign> (*this, mc);
571             return *this;
572         }
573         template<class AE>
574         BOOST_UBLAS_INLINE
operator =(const vector_expression<AE> & ae)575         matrix_column &operator = (const vector_expression<AE> &ae) {
576             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (ae));
577             return *this;
578         }
579         template<class AE>
580         BOOST_UBLAS_INLINE
assign(const vector_expression<AE> & ae)581         matrix_column &assign (const vector_expression<AE> &ae) {
582             vector_assign<scalar_assign> (*this, ae);
583             return *this;
584         }
585         template<class AE>
586         BOOST_UBLAS_INLINE
operator +=(const vector_expression<AE> & ae)587         matrix_column &operator += (const vector_expression<AE> &ae) {
588             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this + ae));
589             return *this;
590         }
591         template<class AE>
592         BOOST_UBLAS_INLINE
plus_assign(const vector_expression<AE> & ae)593         matrix_column &plus_assign (const vector_expression<AE> &ae) {
594             vector_assign<scalar_plus_assign> (*this, ae);
595             return *this;
596         }
597         template<class AE>
598         BOOST_UBLAS_INLINE
operator -=(const vector_expression<AE> & ae)599         matrix_column &operator -= (const vector_expression<AE> &ae) {
600             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this - ae));
601             return *this;
602         }
603         template<class AE>
604         BOOST_UBLAS_INLINE
minus_assign(const vector_expression<AE> & ae)605         matrix_column &minus_assign (const vector_expression<AE> &ae) {
606             vector_assign<scalar_minus_assign> (*this, ae);
607             return *this;
608         }
609         template<class AT>
610         BOOST_UBLAS_INLINE
operator *=(const AT & at)611         matrix_column &operator *= (const AT &at) {
612             vector_assign_scalar<scalar_multiplies_assign> (*this, at);
613             return *this;
614         }
615         template<class AT>
616         BOOST_UBLAS_INLINE
operator /=(const AT & at)617         matrix_column &operator /= (const AT &at) {
618             vector_assign_scalar<scalar_divides_assign> (*this, at);
619             return *this;
620         }
621 
622         // Closure comparison
623         BOOST_UBLAS_INLINE
same_closure(const matrix_column & mc) const624         bool same_closure (const matrix_column &mc) const {
625             return (*this).data_.same_closure (mc.data_);
626         }
627 
628         // Comparison
629         BOOST_UBLAS_INLINE
operator ==(const matrix_column & mc) const630         bool operator == (const matrix_column &mc) const {
631             return (*this).data_ == mc.data_ && index () == mc.index ();
632         }
633 
634         // Swapping
635         BOOST_UBLAS_INLINE
swap(matrix_column mc)636         void swap (matrix_column mc) {
637             if (this != &mc) {
638                 BOOST_UBLAS_CHECK (size () == mc.size (), bad_size ());
639                 // Sparse ranges may be nonconformant now.
640                 // std::swap_ranges (begin (), end (), mc.begin ());
641                 vector_swap<scalar_swap> (*this, mc);
642             }
643         }
644         BOOST_UBLAS_INLINE
swap(matrix_column mc1,matrix_column mc2)645         friend void swap (matrix_column mc1, matrix_column mc2) {
646             mc1.swap (mc2);
647         }
648 
649         // Iterator types
650     private:
651         typedef typename M::const_iterator1 const_subiterator_type;
652         typedef typename boost::mpl::if_<boost::is_const<M>,
653                                           typename M::const_iterator1,
654                                           typename M::iterator1>::type subiterator_type;
655 
656     public:
657 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
658         typedef indexed_iterator<matrix_column<matrix_type>,
659                                  typename subiterator_type::iterator_category> iterator;
660         typedef indexed_const_iterator<matrix_column<matrix_type>,
661                                        typename const_subiterator_type::iterator_category> const_iterator;
662 #else
663         class const_iterator;
664         class iterator;
665 #endif
666 
667         // Element lookup
668         BOOST_UBLAS_INLINE
find(size_type i) const669         const_iterator find (size_type i) const {
670             const_subiterator_type it1 (data_.find1 (1, i, j_));
671 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
672             return const_iterator (*this, it1.index1 ());
673 #else
674             return const_iterator (*this, it1);
675 #endif
676         }
677         BOOST_UBLAS_INLINE
find(size_type i)678         iterator find (size_type i) {
679             subiterator_type it1 (data_.find1 (1, i, j_));
680 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
681             return iterator (*this, it1.index1 ());
682 #else
683             return iterator (*this, it1);
684 #endif
685         }
686 
687 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
688         class const_iterator:
689             public container_const_reference<matrix_column>,
690             public iterator_base_traits<typename const_subiterator_type::iterator_category>::template
691                         iterator_base<const_iterator, value_type>::type {
692         public:
693             typedef typename const_subiterator_type::value_type value_type;
694             typedef typename const_subiterator_type::difference_type difference_type;
695             typedef typename const_subiterator_type::reference reference;
696             typedef typename const_subiterator_type::pointer pointer;
697 
698             // Construction and destruction
699             BOOST_UBLAS_INLINE
const_iterator()700             const_iterator ():
701                 container_const_reference<self_type> (), it_ () {}
702             BOOST_UBLAS_INLINE
const_iterator(const self_type & mc,const const_subiterator_type & it)703             const_iterator (const self_type &mc, const const_subiterator_type &it):
704                 container_const_reference<self_type> (mc), it_ (it) {}
705             BOOST_UBLAS_INLINE
const_iterator(const iterator & it)706             const_iterator (const iterator &it):
707                 container_const_reference<self_type> (it ()), it_ (it.it_) {}
708 
709             // Arithmetic
710             BOOST_UBLAS_INLINE
operator ++()711             const_iterator &operator ++ () {
712                 ++ it_;
713                 return *this;
714             }
715             BOOST_UBLAS_INLINE
operator --()716             const_iterator &operator -- () {
717                 -- it_;
718                 return *this;
719             }
720             BOOST_UBLAS_INLINE
operator +=(difference_type n)721             const_iterator &operator += (difference_type n) {
722                 it_ += n;
723                 return *this;
724             }
725             BOOST_UBLAS_INLINE
operator -=(difference_type n)726             const_iterator &operator -= (difference_type n) {
727                 it_ -= n;
728                 return *this;
729             }
730             BOOST_UBLAS_INLINE
operator -(const const_iterator & it) const731             difference_type operator - (const const_iterator &it) const {
732                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
733                 return it_ - it.it_;
734             }
735 
736             // Dereference
737             BOOST_UBLAS_INLINE
operator *() const738             const_reference operator * () const {
739                 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ());
740                 return *it_;
741             }
742 
743             // Index
744             BOOST_UBLAS_INLINE
index() const745             size_type index () const {
746                 return it_.index1 ();
747             }
748 
749             // Assignment
750             BOOST_UBLAS_INLINE
operator =(const const_iterator & it)751             const_iterator &operator = (const const_iterator &it) {
752                 container_const_reference<self_type>::assign (&it ());
753                 it_ = it.it_;
754                 return *this;
755             }
756 
757             // Comparison
758             BOOST_UBLAS_INLINE
operator ==(const const_iterator & it) const759             bool operator == (const const_iterator &it) const {
760                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
761                 return it_ == it.it_;
762             }
763             BOOST_UBLAS_INLINE
operator <(const const_iterator & it) const764             bool operator < (const const_iterator &it) const {
765                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
766                 return it_ < it.it_;
767             }
768 
769         private:
770             const_subiterator_type it_;
771         };
772 #endif
773 
774         BOOST_UBLAS_INLINE
begin() const775         const_iterator begin () const {
776             return find (0);
777         }
778         BOOST_UBLAS_INLINE
end() const779         const_iterator end () const {
780             return find (size ());
781         }
782 
783 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
784         class iterator:
785             public container_reference<matrix_column>,
786             public iterator_base_traits<typename subiterator_type::iterator_category>::template
787                         iterator_base<iterator, value_type>::type {
788         public:
789             typedef typename subiterator_type::value_type value_type;
790             typedef typename subiterator_type::difference_type difference_type;
791             typedef typename subiterator_type::reference reference;
792             typedef typename subiterator_type::pointer pointer;
793 
794             // Construction and destruction
795             BOOST_UBLAS_INLINE
iterator()796             iterator ():
797                 container_reference<self_type> (), it_ () {}
798             BOOST_UBLAS_INLINE
iterator(self_type & mc,const subiterator_type & it)799             iterator (self_type &mc, const subiterator_type &it):
800                 container_reference<self_type> (mc), it_ (it) {}
801 
802             // Arithmetic
803             BOOST_UBLAS_INLINE
operator ++()804             iterator &operator ++ () {
805                 ++ it_;
806                 return *this;
807             }
808             BOOST_UBLAS_INLINE
operator --()809             iterator &operator -- () {
810                 -- it_;
811                 return *this;
812             }
813             BOOST_UBLAS_INLINE
operator +=(difference_type n)814             iterator &operator += (difference_type n) {
815                 it_ += n;
816                 return *this;
817             }
818             BOOST_UBLAS_INLINE
operator -=(difference_type n)819             iterator &operator -= (difference_type n) {
820                 it_ -= n;
821                 return *this;
822             }
823             BOOST_UBLAS_INLINE
operator -(const iterator & it) const824             difference_type operator - (const iterator &it) const {
825                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
826                 return it_ - it.it_;
827             }
828 
829             // Dereference
830             BOOST_UBLAS_INLINE
operator *() const831             reference operator * () const {
832                 BOOST_UBLAS_CHECK (index () < (*this) ().size (), bad_index ());
833                 return *it_;
834             }
835 
836             // Index
837             BOOST_UBLAS_INLINE
index() const838             size_type index () const {
839                 return it_.index1 ();
840             }
841 
842             // Assignment
843             BOOST_UBLAS_INLINE
operator =(const iterator & it)844             iterator &operator = (const iterator &it) {
845                 container_reference<self_type>::assign (&it ());
846                 it_ = it.it_;
847                 return *this;
848             }
849 
850             // Comparison
851             BOOST_UBLAS_INLINE
operator ==(const iterator & it) const852             bool operator == (const iterator &it) const {
853                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
854                 return it_ == it.it_;
855             }
856             BOOST_UBLAS_INLINE
operator <(const iterator & it) const857             bool operator < (const iterator &it) const {
858                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
859                 return it_ < it.it_;
860             }
861 
862         private:
863             subiterator_type it_;
864 
865             friend class const_iterator;
866         };
867 #endif
868 
869         BOOST_UBLAS_INLINE
begin()870         iterator begin () {
871             return find (0);
872         }
873         BOOST_UBLAS_INLINE
end()874         iterator end () {
875             return find (size ());
876         }
877 
878         // Reverse iterator
879         typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
880         typedef reverse_iterator_base<iterator> reverse_iterator;
881 
882         BOOST_UBLAS_INLINE
rbegin() const883         const_reverse_iterator rbegin () const {
884             return const_reverse_iterator (end ());
885         }
886         BOOST_UBLAS_INLINE
rend() const887         const_reverse_iterator rend () const {
888             return const_reverse_iterator (begin ());
889         }
rbegin()890         reverse_iterator rbegin () {
891             return reverse_iterator (end ());
892         }
893         BOOST_UBLAS_INLINE
rend()894         reverse_iterator rend () {
895             return reverse_iterator (begin ());
896         }
897 
898     private:
899         matrix_closure_type data_;
900         size_type j_;
901     };
902 
903     // Projections
904     template<class M>
905     BOOST_UBLAS_INLINE
column(M & data,typename M::size_type j)906     matrix_column<M> column (M &data, typename M::size_type j) {
907         return matrix_column<M> (data, j);
908     }
909     template<class M>
910     BOOST_UBLAS_INLINE
column(const M & data,typename M::size_type j)911     const matrix_column<const M> column (const M &data, typename M::size_type j) {
912         return matrix_column<const M> (data, j);
913     }
914 
915     // Specialize temporary
916     template <class M>
917     struct vector_temporary_traits< matrix_column<M> >
918     : vector_temporary_traits< M > {} ;
919 
920     // Matrix based vector range class
921     template<class M>
922     class matrix_vector_range:
923         public vector_expression<matrix_vector_range<M> > {
924 
925         typedef matrix_vector_range<M> self_type;
926     public:
927 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
928         using vector_expression<self_type>::operator ();
929 #endif
930         typedef M matrix_type;
931         typedef typename M::size_type size_type;
932         typedef typename M::difference_type difference_type;
933         typedef typename M::value_type value_type;
934         typedef typename M::const_reference const_reference;
935         typedef typename boost::mpl::if_<boost::is_const<M>,
936                                           typename M::const_reference,
937                                           typename M::reference>::type reference;
938         typedef typename boost::mpl::if_<boost::is_const<M>,
939                                           typename M::const_closure_type,
940                                           typename M::closure_type>::type matrix_closure_type;
941         typedef basic_range<size_type, difference_type> range_type;
942         typedef const self_type const_closure_type;
943         typedef self_type closure_type;
944         typedef typename storage_restrict_traits<typename M::storage_category,
945                                                  dense_proxy_tag>::storage_category storage_category;
946 
947         // Construction and destruction
948         BOOST_UBLAS_INLINE
matrix_vector_range(matrix_type & data,const range_type & r1,const range_type & r2)949         matrix_vector_range (matrix_type &data, const range_type &r1, const range_type &r2):
950             data_ (data), r1_ (r1.preprocess (data.size1 ())), r2_ (r2.preprocess (data.size2 ())) {
951             // Early checking of preconditions here.
952             // BOOST_UBLAS_CHECK (r1_.start () <= data_.size1 () &&
953             //                     r1_.start () + r1_.size () <= data_.size1 (), bad_index ());
954             // BOOST_UBLAS_CHECK (r2_.start () <= data_.size2 () &&
955             //                   r2_.start () + r2_.size () <= data_.size2 (), bad_index ());
956             // BOOST_UBLAS_CHECK (r1_.size () == r2_.size (), bad_size ());
957         }
958 
959         // Accessors
960         BOOST_UBLAS_INLINE
start1() const961         size_type start1 () const {
962             return r1_.start ();
963         }
964         BOOST_UBLAS_INLINE
start2() const965         size_type start2 () const {
966             return r2_.start ();
967         }
968         BOOST_UBLAS_INLINE
size() const969         size_type size () const {
970             return BOOST_UBLAS_SAME (r1_.size (), r2_.size ());
971         }
972 
973         // Storage accessors
974         BOOST_UBLAS_INLINE
data() const975         const matrix_closure_type &data () const {
976             return data_;
977         }
978         BOOST_UBLAS_INLINE
data()979         matrix_closure_type &data () {
980             return data_;
981         }
982 
983         // Element access
984 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
985         BOOST_UBLAS_INLINE
operator ()(size_type i) const986         const_reference operator () (size_type i) const {
987             return data_ (r1_ (i), r2_ (i));
988         }
989         BOOST_UBLAS_INLINE
operator ()(size_type i)990         reference operator () (size_type i) {
991             return data_ (r1_ (i), r2_ (i));
992         }
993 
994         BOOST_UBLAS_INLINE
operator [](size_type i) const995         const_reference operator [] (size_type i) const {
996             return (*this) (i);
997         }
998         BOOST_UBLAS_INLINE
operator [](size_type i)999         reference operator [] (size_type i) {
1000             return (*this) (i);
1001         }
1002 #else
1003         BOOST_UBLAS_INLINE
operator ()(size_type i) const1004         reference operator () (size_type i) const {
1005             return data_ (r1_ (i), r2_ (i));
1006         }
1007 
1008         BOOST_UBLAS_INLINE
operator [](size_type i) const1009         reference operator [] (size_type i) const {
1010             return (*this) (i);
1011         }
1012 #endif
1013 
1014         // Assignment
1015         BOOST_UBLAS_INLINE
operator =(const matrix_vector_range & mvr)1016         matrix_vector_range &operator = (const matrix_vector_range &mvr) {
1017             // ISSUE need a temporary, proxy can be overlaping alias
1018             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (mvr));
1019             return *this;
1020         }
1021         BOOST_UBLAS_INLINE
assign_temporary(matrix_vector_range & mvr)1022         matrix_vector_range &assign_temporary (matrix_vector_range &mvr) {
1023             // assign elements, proxied container remains the same
1024             vector_assign<scalar_assign> (*this, mvr);
1025             return *this;
1026         }
1027         template<class AE>
1028         BOOST_UBLAS_INLINE
operator =(const vector_expression<AE> & ae)1029         matrix_vector_range &operator = (const vector_expression<AE> &ae) {
1030             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (ae));
1031             return *this;
1032         }
1033         template<class AE>
1034         BOOST_UBLAS_INLINE
assign(const vector_expression<AE> & ae)1035         matrix_vector_range &assign (const vector_expression<AE> &ae) {
1036             vector_assign<scalar_assign> (*this, ae);
1037             return *this;
1038         }
1039         template<class AE>
1040         BOOST_UBLAS_INLINE
operator +=(const vector_expression<AE> & ae)1041         matrix_vector_range &operator += (const vector_expression<AE> &ae) {
1042             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this + ae));
1043             return *this;
1044         }
1045         template<class AE>
1046         BOOST_UBLAS_INLINE
plus_assign(const vector_expression<AE> & ae)1047         matrix_vector_range &plus_assign (const vector_expression<AE> &ae) {
1048             vector_assign<scalar_plus_assign> (*this, ae);
1049             return *this;
1050         }
1051         template<class AE>
1052         BOOST_UBLAS_INLINE
operator -=(const vector_expression<AE> & ae)1053         matrix_vector_range &operator -= (const vector_expression<AE> &ae) {
1054             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this - ae));
1055             return *this;
1056         }
1057         template<class AE>
1058         BOOST_UBLAS_INLINE
minus_assign(const vector_expression<AE> & ae)1059         matrix_vector_range &minus_assign (const vector_expression<AE> &ae) {
1060             vector_assign<scalar_minus_assign> (*this, ae);
1061             return *this;
1062         }
1063         template<class AT>
1064         BOOST_UBLAS_INLINE
operator *=(const AT & at)1065         matrix_vector_range &operator *= (const AT &at) {
1066             vector_assign_scalar<scalar_multiplies_assign> (*this, at);
1067             return *this;
1068         }
1069         template<class AT>
1070         BOOST_UBLAS_INLINE
operator /=(const AT & at)1071         matrix_vector_range &operator /= (const AT &at) {
1072             vector_assign_scalar<scalar_divides_assign> (*this, at);
1073             return *this;
1074         }
1075 
1076         // Closure comparison
1077         BOOST_UBLAS_INLINE
same_closure(const matrix_vector_range & mvr) const1078         bool same_closure (const matrix_vector_range &mvr) const {
1079             return (*this).data_.same_closure (mvr.data_);
1080         }
1081 
1082         // Comparison
1083         BOOST_UBLAS_INLINE
operator ==(const matrix_vector_range & mvr) const1084         bool operator == (const matrix_vector_range &mvr) const {
1085             return (*this).data_ == mvr.data_ && r1_ == mvr.r1_ && r2_ == mvr.r2_;
1086         }
1087 
1088         // Swapping
1089         BOOST_UBLAS_INLINE
swap(matrix_vector_range mvr)1090         void swap (matrix_vector_range mvr) {
1091             if (this != &mvr) {
1092                 BOOST_UBLAS_CHECK (size () == mvr.size (), bad_size ());
1093                 // Sparse ranges may be nonconformant now.
1094                 // std::swap_ranges (begin (), end (), mvr.begin ());
1095                 vector_swap<scalar_swap> (*this, mvr);
1096             }
1097         }
1098         BOOST_UBLAS_INLINE
swap(matrix_vector_range mvr1,matrix_vector_range mvr2)1099         friend void swap (matrix_vector_range mvr1, matrix_vector_range mvr2) {
1100             mvr1.swap (mvr2);
1101         }
1102 
1103         // Iterator types
1104     private:
1105         // Use range as an index - FIXME this fails for packed assignment
1106         typedef typename range_type::const_iterator const_subiterator1_type;
1107         typedef typename range_type::const_iterator subiterator1_type;
1108         typedef typename range_type::const_iterator const_subiterator2_type;
1109         typedef typename range_type::const_iterator subiterator2_type;
1110 
1111     public:
1112         class const_iterator;
1113         class iterator;
1114 
1115         // Element lookup
1116         BOOST_UBLAS_INLINE
find(size_type i) const1117         const_iterator find (size_type i) const {
1118             return const_iterator (*this, r1_.begin () + i, r2_.begin () + i);
1119         }
1120         BOOST_UBLAS_INLINE
find(size_type i)1121         iterator find (size_type i) {
1122             return iterator (*this, r1_.begin () + i, r2_.begin () + i);
1123         }
1124 
1125         class const_iterator:
1126             public container_const_reference<matrix_vector_range>,
1127             public iterator_base_traits<typename M::const_iterator1::iterator_category>::template
1128                         iterator_base<const_iterator, value_type>::type {
1129         public:
1130             // FIXME Iterator can never be different code was:
1131             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
1132             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
1133 
1134             typedef typename matrix_vector_range::value_type value_type;
1135             typedef typename matrix_vector_range::difference_type difference_type;
1136             typedef typename matrix_vector_range::const_reference reference;
1137             typedef const typename matrix_vector_range::value_type *pointer;
1138 
1139             // Construction and destruction
1140             BOOST_UBLAS_INLINE
const_iterator()1141             const_iterator ():
1142                 container_const_reference<self_type> (), it1_ (), it2_ () {}
1143             BOOST_UBLAS_INLINE
const_iterator(const self_type & mvr,const const_subiterator1_type & it1,const const_subiterator2_type & it2)1144             const_iterator (const self_type &mvr, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
1145                 container_const_reference<self_type> (mvr), it1_ (it1), it2_ (it2) {}
1146             BOOST_UBLAS_INLINE
const_iterator(const iterator & it)1147             const_iterator (const iterator &it):
1148                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
1149 
1150             // Arithmetic
1151             BOOST_UBLAS_INLINE
operator ++()1152             const_iterator &operator ++ () {
1153                 ++ it1_;
1154                 ++ it2_;
1155                 return *this;
1156             }
1157             BOOST_UBLAS_INLINE
operator --()1158             const_iterator &operator -- () {
1159                 -- it1_;
1160                 -- it2_;
1161                 return *this;
1162             }
1163             BOOST_UBLAS_INLINE
operator +=(difference_type n)1164             const_iterator &operator += (difference_type n) {
1165                 it1_ += n;
1166                 it2_ += n;
1167                 return *this;
1168             }
1169             BOOST_UBLAS_INLINE
operator -=(difference_type n)1170             const_iterator &operator -= (difference_type n) {
1171                 it1_ -= n;
1172                 it2_ -= n;
1173                 return *this;
1174             }
1175             BOOST_UBLAS_INLINE
operator -(const const_iterator & it) const1176             difference_type operator - (const const_iterator &it) const {
1177                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
1178                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
1179             }
1180 
1181             // Dereference
1182             BOOST_UBLAS_INLINE
operator *() const1183             const_reference operator * () const {
1184                 // FIXME replace find with at_element
1185                 return (*this) ().data_ (*it1_, *it2_);
1186             }
1187 
1188             // Index
1189             BOOST_UBLAS_INLINE
index() const1190             size_type  index () const {
1191                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
1192             }
1193 
1194             // Assignment
1195             BOOST_UBLAS_INLINE
operator =(const const_iterator & it)1196             const_iterator &operator = (const const_iterator &it) {
1197                 container_const_reference<self_type>::assign (&it ());
1198                 it1_ = it.it1_;
1199                 it2_ = it.it2_;
1200                 return *this;
1201             }
1202 
1203             // Comparison
1204             BOOST_UBLAS_INLINE
operator ==(const const_iterator & it) const1205             bool operator == (const const_iterator &it) const {
1206                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1207                 return it1_ == it.it1_ && it2_ == it.it2_;
1208             }
1209             BOOST_UBLAS_INLINE
operator <(const const_iterator & it) const1210             bool operator < (const const_iterator &it) const {
1211                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1212                 return it1_ < it.it1_ && it2_ < it.it2_;
1213             }
1214 
1215         private:
1216             const_subiterator1_type it1_;
1217             const_subiterator2_type it2_;
1218         };
1219 
1220         BOOST_UBLAS_INLINE
begin() const1221         const_iterator begin () const {
1222             return find (0);
1223         }
1224         BOOST_UBLAS_INLINE
end() const1225         const_iterator end () const {
1226             return find (size ());
1227         }
1228 
1229         class iterator:
1230             public container_reference<matrix_vector_range>,
1231             public iterator_base_traits<typename M::iterator1::iterator_category>::template
1232                         iterator_base<iterator, value_type>::type {
1233         public:
1234             // FIXME Iterator can never be different code was:
1235             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
1236             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
1237 
1238             typedef typename matrix_vector_range::value_type value_type;
1239             typedef typename matrix_vector_range::difference_type difference_type;
1240             typedef typename matrix_vector_range::reference reference;
1241             typedef typename matrix_vector_range::value_type *pointer;
1242 
1243             // Construction and destruction
1244             BOOST_UBLAS_INLINE
iterator()1245             iterator ():
1246                 container_reference<self_type> (), it1_ (), it2_ () {}
1247             BOOST_UBLAS_INLINE
iterator(self_type & mvr,const subiterator1_type & it1,const subiterator2_type & it2)1248             iterator (self_type &mvr, const subiterator1_type &it1, const subiterator2_type &it2):
1249                 container_reference<self_type> (mvr), it1_ (it1), it2_ (it2) {}
1250 
1251             // Arithmetic
1252             BOOST_UBLAS_INLINE
operator ++()1253             iterator &operator ++ () {
1254                 ++ it1_;
1255                 ++ it2_;
1256                 return *this;
1257             }
1258             BOOST_UBLAS_INLINE
operator --()1259             iterator &operator -- () {
1260                 -- it1_;
1261                 -- it2_;
1262                 return *this;
1263             }
1264             BOOST_UBLAS_INLINE
operator +=(difference_type n)1265             iterator &operator += (difference_type n) {
1266                 it1_ += n;
1267                 it2_ += n;
1268                 return *this;
1269             }
1270             BOOST_UBLAS_INLINE
operator -=(difference_type n)1271             iterator &operator -= (difference_type n) {
1272                 it1_ -= n;
1273                 it2_ -= n;
1274                 return *this;
1275             }
1276             BOOST_UBLAS_INLINE
operator -(const iterator & it) const1277             difference_type operator - (const iterator &it) const {
1278                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
1279                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
1280             }
1281 
1282             // Dereference
1283             BOOST_UBLAS_INLINE
operator *() const1284             reference operator * () const {
1285                 // FIXME replace find with at_element
1286                 return (*this) ().data_ (*it1_, *it2_);
1287             }
1288 
1289             // Index
1290             BOOST_UBLAS_INLINE
index() const1291             size_type index () const {
1292                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
1293             }
1294 
1295             // Assignment
1296             BOOST_UBLAS_INLINE
operator =(const iterator & it)1297             iterator &operator = (const iterator &it) {
1298                 container_reference<self_type>::assign (&it ());
1299                 it1_ = it.it1_;
1300                 it2_ = it.it2_;
1301                 return *this;
1302             }
1303 
1304             // Comparison
1305             BOOST_UBLAS_INLINE
operator ==(const iterator & it) const1306             bool operator == (const iterator &it) const {
1307                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1308                 return it1_ == it.it1_ && it2_ == it.it2_;
1309             }
1310             BOOST_UBLAS_INLINE
operator <(const iterator & it) const1311             bool operator < (const iterator &it) const {
1312                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1313                 return it1_ < it.it1_ && it2_ < it.it2_;
1314             }
1315 
1316         private:
1317             subiterator1_type it1_;
1318             subiterator2_type it2_;
1319 
1320             friend class const_iterator;
1321         };
1322 
1323         BOOST_UBLAS_INLINE
begin()1324         iterator begin () {
1325             return find (0);
1326         }
1327         BOOST_UBLAS_INLINE
end()1328         iterator end () {
1329             return find (size ());
1330         }
1331 
1332         // Reverse iterator
1333         typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
1334         typedef reverse_iterator_base<iterator> reverse_iterator;
1335 
1336         BOOST_UBLAS_INLINE
rbegin() const1337         const_reverse_iterator rbegin () const {
1338             return const_reverse_iterator (end ());
1339         }
1340         BOOST_UBLAS_INLINE
rend() const1341         const_reverse_iterator rend () const {
1342             return const_reverse_iterator (begin ());
1343         }
1344         BOOST_UBLAS_INLINE
rbegin()1345         reverse_iterator rbegin () {
1346             return reverse_iterator (end ());
1347         }
1348         BOOST_UBLAS_INLINE
rend()1349         reverse_iterator rend () {
1350             return reverse_iterator (begin ());
1351         }
1352 
1353     private:
1354         matrix_closure_type data_;
1355         range_type r1_;
1356         range_type r2_;
1357     };
1358 
1359     // Specialize temporary
1360     template <class M>
1361     struct vector_temporary_traits< matrix_vector_range<M> >
1362     : vector_temporary_traits< M > {} ;
1363 
1364     // Matrix based vector slice class
1365     template<class M>
1366     class matrix_vector_slice:
1367         public vector_expression<matrix_vector_slice<M> > {
1368 
1369         typedef matrix_vector_slice<M> self_type;
1370     public:
1371 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1372         using vector_expression<self_type>::operator ();
1373 #endif
1374         typedef M matrix_type;
1375         typedef typename M::size_type size_type;
1376         typedef typename M::difference_type difference_type;
1377         typedef typename M::value_type value_type;
1378         typedef typename M::const_reference const_reference;
1379         typedef typename boost::mpl::if_<boost::is_const<M>,
1380                                           typename M::const_reference,
1381                                           typename M::reference>::type reference;
1382         typedef typename boost::mpl::if_<boost::is_const<M>,
1383                                           typename M::const_closure_type,
1384                                           typename M::closure_type>::type matrix_closure_type;
1385         typedef basic_range<size_type, difference_type> range_type;
1386         typedef basic_slice<size_type, difference_type> slice_type;
1387         typedef const self_type const_closure_type;
1388         typedef self_type closure_type;
1389         typedef typename storage_restrict_traits<typename M::storage_category,
1390                                                  dense_proxy_tag>::storage_category storage_category;
1391 
1392         // Construction and destruction
1393         BOOST_UBLAS_INLINE
matrix_vector_slice(matrix_type & data,const slice_type & s1,const slice_type & s2)1394         matrix_vector_slice (matrix_type &data, const slice_type &s1, const slice_type &s2):
1395             data_ (data), s1_ (s1.preprocess (data.size1 ())), s2_ (s2.preprocess (data.size2 ())) {
1396             // Early checking of preconditions here.
1397             // BOOST_UBLAS_CHECK (s1_.start () <= data_.size1 () &&
1398             //                    s1_.start () + s1_.stride () * (s1_.size () - (s1_.size () > 0)) <= data_.size1 (), bad_index ());
1399             // BOOST_UBLAS_CHECK (s2_.start () <= data_.size2 () &&
1400             //                    s2_.start () + s2_.stride () * (s2_.size () - (s2_.size () > 0)) <= data_.size2 (), bad_index ());
1401         }
1402 
1403         // Accessors
1404         BOOST_UBLAS_INLINE
start1() const1405         size_type start1 () const {
1406             return s1_.start ();
1407         }
1408         BOOST_UBLAS_INLINE
start2() const1409         size_type start2 () const {
1410             return s2_.start ();
1411         }
1412         BOOST_UBLAS_INLINE
stride1() const1413         difference_type stride1 () const {
1414             return s1_.stride ();
1415         }
1416         BOOST_UBLAS_INLINE
stride2() const1417         difference_type stride2 () const {
1418             return s2_.stride ();
1419         }
1420         BOOST_UBLAS_INLINE
size() const1421         size_type size () const {
1422             return BOOST_UBLAS_SAME (s1_.size (), s2_.size ());
1423         }
1424 
1425         // Storage accessors
1426         BOOST_UBLAS_INLINE
data() const1427         const matrix_closure_type &data () const {
1428             return data_;
1429         }
1430         BOOST_UBLAS_INLINE
data()1431         matrix_closure_type &data () {
1432             return data_;
1433         }
1434 
1435         // Element access
1436 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
1437         BOOST_UBLAS_INLINE
operator ()(size_type i) const1438         const_reference operator () (size_type i) const {
1439             return data_ (s1_ (i), s2_ (i));
1440         }
1441         BOOST_UBLAS_INLINE
operator ()(size_type i)1442         reference operator () (size_type i) {
1443             return data_ (s1_ (i), s2_ (i));
1444         }
1445 
1446         BOOST_UBLAS_INLINE
operator [](size_type i) const1447         const_reference operator [] (size_type i) const {
1448             return (*this) (i);
1449         }
1450         BOOST_UBLAS_INLINE
operator [](size_type i)1451         reference operator [] (size_type i) {
1452             return (*this) (i);
1453         }
1454 #else
1455         BOOST_UBLAS_INLINE
operator ()(size_type i) const1456         reference operator () (size_type i) const {
1457             return data_ (s1_ (i), s2_ (i));
1458         }
1459 
1460         BOOST_UBLAS_INLINE
operator [](size_type i) const1461         reference operator [] (size_type i) const {
1462             return (*this) (i);
1463         }
1464 #endif
1465 
1466         // Assignment
1467         BOOST_UBLAS_INLINE
operator =(const matrix_vector_slice & mvs)1468         matrix_vector_slice &operator = (const matrix_vector_slice &mvs) {
1469             // ISSUE need a temporary, proxy can be overlaping alias
1470             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (mvs));
1471             return *this;
1472         }
1473         BOOST_UBLAS_INLINE
assign_temporary(matrix_vector_slice & mvs)1474         matrix_vector_slice &assign_temporary (matrix_vector_slice &mvs) {
1475             // assign elements, proxied container remains the same
1476             vector_assign<scalar_assign> (*this, mvs);
1477             return *this;
1478         }
1479         template<class AE>
1480         BOOST_UBLAS_INLINE
operator =(const vector_expression<AE> & ae)1481         matrix_vector_slice &operator = (const vector_expression<AE> &ae) {
1482             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (ae));
1483             return *this;
1484         }
1485         template<class AE>
1486         BOOST_UBLAS_INLINE
assign(const vector_expression<AE> & ae)1487         matrix_vector_slice &assign (const vector_expression<AE> &ae) {
1488             vector_assign<scalar_assign> (*this, ae);
1489             return *this;
1490         }
1491         template<class AE>
1492         BOOST_UBLAS_INLINE
operator +=(const vector_expression<AE> & ae)1493         matrix_vector_slice &operator += (const vector_expression<AE> &ae) {
1494             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this + ae));
1495             return *this;
1496         }
1497         template<class AE>
1498         BOOST_UBLAS_INLINE
plus_assign(const vector_expression<AE> & ae)1499         matrix_vector_slice &plus_assign (const vector_expression<AE> &ae) {
1500             vector_assign<scalar_plus_assign> (*this, ae);
1501             return *this;
1502         }
1503         template<class AE>
1504         BOOST_UBLAS_INLINE
operator -=(const vector_expression<AE> & ae)1505         matrix_vector_slice &operator -= (const vector_expression<AE> &ae) {
1506             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this - ae));
1507             return *this;
1508         }
1509         template<class AE>
1510         BOOST_UBLAS_INLINE
minus_assign(const vector_expression<AE> & ae)1511         matrix_vector_slice &minus_assign (const vector_expression<AE> &ae) {
1512             vector_assign<scalar_minus_assign> (*this, ae);
1513             return *this;
1514         }
1515         template<class AT>
1516         BOOST_UBLAS_INLINE
operator *=(const AT & at)1517         matrix_vector_slice &operator *= (const AT &at) {
1518             vector_assign_scalar<scalar_multiplies_assign> (*this, at);
1519             return *this;
1520         }
1521         template<class AT>
1522         BOOST_UBLAS_INLINE
operator /=(const AT & at)1523         matrix_vector_slice &operator /= (const AT &at) {
1524             vector_assign_scalar<scalar_divides_assign> (*this, at);
1525             return *this;
1526         }
1527 
1528         // Closure comparison
1529         BOOST_UBLAS_INLINE
same_closure(const matrix_vector_slice & mvs) const1530         bool same_closure (const matrix_vector_slice &mvs) const {
1531             return (*this).data_.same_closure (mvs.data_);
1532         }
1533 
1534         // Comparison
1535         BOOST_UBLAS_INLINE
operator ==(const matrix_vector_slice & mvs) const1536         bool operator == (const matrix_vector_slice &mvs) const {
1537             return (*this).data_ == mvs.data_ && s1_ == mvs.s1_ && s2_ == mvs.s2_;
1538         }
1539 
1540         // Swapping
1541         BOOST_UBLAS_INLINE
swap(matrix_vector_slice mvs)1542         void swap (matrix_vector_slice mvs) {
1543             if (this != &mvs) {
1544                 BOOST_UBLAS_CHECK (size () == mvs.size (), bad_size ());
1545                 // Sparse ranges may be nonconformant now.
1546                 // std::swap_ranges (begin (), end (), mvs.begin ());
1547                 vector_swap<scalar_swap> (*this, mvs);
1548             }
1549         }
1550         BOOST_UBLAS_INLINE
swap(matrix_vector_slice mvs1,matrix_vector_slice mvs2)1551         friend void swap (matrix_vector_slice mvs1, matrix_vector_slice mvs2) {
1552             mvs1.swap (mvs2);
1553         }
1554 
1555         // Iterator types
1556     private:
1557         // Use slice as an index - FIXME this fails for packed assignment
1558         typedef typename slice_type::const_iterator const_subiterator1_type;
1559         typedef typename slice_type::const_iterator subiterator1_type;
1560         typedef typename slice_type::const_iterator const_subiterator2_type;
1561         typedef typename slice_type::const_iterator subiterator2_type;
1562 
1563     public:
1564         class const_iterator;
1565         class iterator;
1566 
1567         // Element lookup
1568         BOOST_UBLAS_INLINE
find(size_type i) const1569         const_iterator find (size_type i) const {
1570             return const_iterator (*this, s1_.begin () + i, s2_.begin () + i);
1571         }
1572         BOOST_UBLAS_INLINE
find(size_type i)1573         iterator find (size_type i) {
1574             return iterator (*this, s1_.begin () + i, s2_.begin () + i);
1575         }
1576 
1577         // Iterators simply are indices.
1578 
1579         class const_iterator:
1580             public container_const_reference<matrix_vector_slice>,
1581             public iterator_base_traits<typename M::const_iterator1::iterator_category>::template
1582                         iterator_base<const_iterator, value_type>::type {
1583         public:
1584             // FIXME Iterator can never be different code was:
1585             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
1586             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
1587 
1588             typedef typename matrix_vector_slice::value_type value_type;
1589             typedef typename matrix_vector_slice::difference_type difference_type;
1590             typedef typename matrix_vector_slice::const_reference reference;
1591             typedef const typename matrix_vector_slice::value_type *pointer;
1592 
1593             // Construction and destruction
1594             BOOST_UBLAS_INLINE
const_iterator()1595             const_iterator ():
1596                 container_const_reference<self_type> (), it1_ (), it2_ () {}
1597             BOOST_UBLAS_INLINE
const_iterator(const self_type & mvs,const const_subiterator1_type & it1,const const_subiterator2_type & it2)1598             const_iterator (const self_type &mvs, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
1599                 container_const_reference<self_type> (mvs), it1_ (it1), it2_ (it2) {}
1600             BOOST_UBLAS_INLINE
const_iterator(const iterator & it)1601             const_iterator (const iterator &it):
1602                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
1603 
1604             // Arithmetic
1605             BOOST_UBLAS_INLINE
operator ++()1606             const_iterator &operator ++ () {
1607                 ++ it1_;
1608                 ++ it2_;
1609                 return *this;
1610             }
1611             BOOST_UBLAS_INLINE
operator --()1612             const_iterator &operator -- () {
1613                 -- it1_;
1614                 -- it2_;
1615                 return *this;
1616             }
1617             BOOST_UBLAS_INLINE
operator +=(difference_type n)1618             const_iterator &operator += (difference_type n) {
1619                 it1_ += n;
1620                 it2_ += n;
1621                 return *this;
1622             }
1623             BOOST_UBLAS_INLINE
operator -=(difference_type n)1624             const_iterator &operator -= (difference_type n) {
1625                 it1_ -= n;
1626                 it2_ -= n;
1627                 return *this;
1628             }
1629             BOOST_UBLAS_INLINE
operator -(const const_iterator & it) const1630             difference_type operator - (const const_iterator &it) const {
1631                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
1632                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
1633             }
1634 
1635             // Dereference
1636             BOOST_UBLAS_INLINE
operator *() const1637             const_reference operator * () const {
1638                 // FIXME replace find with at_element
1639                 return (*this) ().data_ (*it1_, *it2_);
1640             }
1641 
1642             // Index
1643             BOOST_UBLAS_INLINE
index() const1644             size_type  index () const {
1645                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
1646             }
1647 
1648             // Assignment
1649             BOOST_UBLAS_INLINE
operator =(const const_iterator & it)1650             const_iterator &operator = (const const_iterator &it) {
1651                 container_const_reference<self_type>::assign (&it ());
1652                 it1_ = it.it1_;
1653                 it2_ = it.it2_;
1654                 return *this;
1655             }
1656 
1657             // Comparison
1658             BOOST_UBLAS_INLINE
operator ==(const const_iterator & it) const1659             bool operator == (const const_iterator &it) const {
1660                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1661                 return it1_ == it.it1_ && it2_ == it.it2_;
1662             }
1663             BOOST_UBLAS_INLINE
operator <(const const_iterator & it) const1664             bool operator < (const const_iterator &it) const {
1665                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1666                 return it1_ < it.it1_ && it2_ < it.it2_;
1667             }
1668 
1669         private:
1670             const_subiterator1_type it1_;
1671             const_subiterator2_type it2_;
1672         };
1673 
1674         BOOST_UBLAS_INLINE
begin() const1675         const_iterator begin () const {
1676             return find (0);
1677         }
1678         BOOST_UBLAS_INLINE
end() const1679         const_iterator end () const {
1680             return find (size ());
1681         }
1682 
1683         class iterator:
1684             public container_reference<matrix_vector_slice>,
1685             public iterator_base_traits<typename M::iterator1::iterator_category>::template
1686                         iterator_base<iterator, value_type>::type {
1687         public:
1688             // FIXME Iterator can never be different code was:
1689             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
1690             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
1691 
1692             typedef typename matrix_vector_slice::value_type value_type;
1693             typedef typename matrix_vector_slice::difference_type difference_type;
1694             typedef typename matrix_vector_slice::reference reference;
1695             typedef typename matrix_vector_slice::value_type *pointer;
1696 
1697             // Construction and destruction
1698             BOOST_UBLAS_INLINE
iterator()1699             iterator ():
1700                 container_reference<self_type> (), it1_ (), it2_ () {}
1701             BOOST_UBLAS_INLINE
iterator(self_type & mvs,const subiterator1_type & it1,const subiterator2_type & it2)1702             iterator (self_type &mvs, const subiterator1_type &it1, const subiterator2_type &it2):
1703                 container_reference<self_type> (mvs), it1_ (it1), it2_ (it2) {}
1704 
1705             // Arithmetic
1706             BOOST_UBLAS_INLINE
operator ++()1707             iterator &operator ++ () {
1708                 ++ it1_;
1709                 ++ it2_;
1710                 return *this;
1711             }
1712             BOOST_UBLAS_INLINE
operator --()1713             iterator &operator -- () {
1714                 -- it1_;
1715                 -- it2_;
1716                 return *this;
1717             }
1718             BOOST_UBLAS_INLINE
operator +=(difference_type n)1719             iterator &operator += (difference_type n) {
1720                 it1_ += n;
1721                 it2_ += n;
1722                 return *this;
1723             }
1724             BOOST_UBLAS_INLINE
operator -=(difference_type n)1725             iterator &operator -= (difference_type n) {
1726                 it1_ -= n;
1727                 it2_ -= n;
1728                 return *this;
1729             }
1730             BOOST_UBLAS_INLINE
operator -(const iterator & it) const1731             difference_type operator - (const iterator &it) const {
1732                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
1733                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
1734             }
1735 
1736             // Dereference
1737             BOOST_UBLAS_INLINE
operator *() const1738             reference operator * () const {
1739                 // FIXME replace find with at_element
1740                 return (*this) ().data_ (*it1_, *it2_);
1741             }
1742 
1743             // Index
1744             BOOST_UBLAS_INLINE
index() const1745             size_type index () const {
1746                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
1747             }
1748 
1749             // Assignment
1750             BOOST_UBLAS_INLINE
operator =(const iterator & it)1751             iterator &operator = (const iterator &it) {
1752                 container_reference<self_type>::assign (&it ());
1753                 it1_ = it.it1_;
1754                 it2_ = it.it2_;
1755                 return *this;
1756             }
1757 
1758             // Comparison
1759             BOOST_UBLAS_INLINE
operator ==(const iterator & it) const1760             bool operator == (const iterator &it) const {
1761                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1762                 return it1_ == it.it1_ && it2_ == it.it2_;
1763             }
1764             BOOST_UBLAS_INLINE
operator <(const iterator & it) const1765             bool operator < (const iterator &it) const {
1766                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1767                 return it1_ < it.it1_ && it2_ < it.it2_;
1768             }
1769 
1770         private:
1771             subiterator1_type it1_;
1772             subiterator2_type it2_;
1773 
1774             friend class const_iterator;
1775         };
1776 
1777         BOOST_UBLAS_INLINE
begin()1778         iterator begin () {
1779             return find (0);
1780         }
1781         BOOST_UBLAS_INLINE
end()1782         iterator end () {
1783             return find (size ());
1784         }
1785 
1786         // Reverse iterator
1787         typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
1788         typedef reverse_iterator_base<iterator> reverse_iterator;
1789 
1790         BOOST_UBLAS_INLINE
rbegin() const1791         const_reverse_iterator rbegin () const {
1792             return const_reverse_iterator (end ());
1793         }
1794         BOOST_UBLAS_INLINE
rend() const1795         const_reverse_iterator rend () const {
1796             return const_reverse_iterator (begin ());
1797         }
1798         BOOST_UBLAS_INLINE
rbegin()1799         reverse_iterator rbegin () {
1800             return reverse_iterator (end ());
1801         }
1802         BOOST_UBLAS_INLINE
rend()1803         reverse_iterator rend () {
1804             return reverse_iterator (begin ());
1805         }
1806 
1807     private:
1808         matrix_closure_type data_;
1809         slice_type s1_;
1810         slice_type s2_;
1811     };
1812 
1813     // Specialize temporary
1814     template <class M>
1815     struct vector_temporary_traits< matrix_vector_slice<M> >
1816     : vector_temporary_traits< M > {} ;
1817 
1818     // Matrix based vector indirection class
1819     template<class M, class IA>
1820     class matrix_vector_indirect:
1821         public vector_expression<matrix_vector_indirect<M, IA> > {
1822 
1823         typedef matrix_vector_indirect<M, IA> self_type;
1824     public:
1825 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1826         using vector_expression<self_type>::operator ();
1827 #endif
1828         typedef M matrix_type;
1829         typedef IA indirect_array_type;
1830         typedef typename M::size_type size_type;
1831         typedef typename M::difference_type difference_type;
1832         typedef typename M::value_type value_type;
1833         typedef typename M::const_reference const_reference;
1834         typedef typename boost::mpl::if_<boost::is_const<M>,
1835                                           typename M::const_reference,
1836                                           typename M::reference>::type reference;
1837         typedef typename boost::mpl::if_<boost::is_const<M>,
1838                                           typename M::const_closure_type,
1839                                           typename M::closure_type>::type matrix_closure_type;
1840         typedef const self_type const_closure_type;
1841         typedef self_type closure_type;
1842         typedef typename storage_restrict_traits<typename M::storage_category,
1843                                                  dense_proxy_tag>::storage_category storage_category;
1844 
1845         // Construction and destruction
1846         BOOST_UBLAS_INLINE
matrix_vector_indirect(matrix_type & data,size_type size)1847         matrix_vector_indirect (matrix_type &data, size_type size):
1848             data_ (data), ia1_ (size), ia2_ (size) {}
1849         BOOST_UBLAS_INLINE
matrix_vector_indirect(matrix_type & data,const indirect_array_type & ia1,const indirect_array_type & ia2)1850         matrix_vector_indirect (matrix_type &data, const indirect_array_type &ia1, const indirect_array_type &ia2):
1851             data_ (data), ia1_ (ia1), ia2_ (ia2) {
1852             // Early checking of preconditions here.
1853             // BOOST_UBLAS_CHECK (ia1_.size () == ia2_.size (), bad_size ());
1854         }
1855 
1856         // Accessors
1857         BOOST_UBLAS_INLINE
size() const1858         size_type size () const {
1859             return BOOST_UBLAS_SAME (ia1_.size (), ia2_.size ());
1860         }
1861         BOOST_UBLAS_INLINE
indirect1() const1862         const indirect_array_type &indirect1 () const {
1863             return ia1_;
1864         }
1865         BOOST_UBLAS_INLINE
indirect1()1866         indirect_array_type &indirect1 () {
1867             return ia1_;
1868         }
1869         BOOST_UBLAS_INLINE
indirect2() const1870         const indirect_array_type &indirect2 () const {
1871             return ia2_;
1872         }
1873         BOOST_UBLAS_INLINE
indirect2()1874         indirect_array_type &indirect2 () {
1875             return ia2_;
1876         }
1877 
1878         // Storage accessors
1879         BOOST_UBLAS_INLINE
data() const1880         const matrix_closure_type &data () const {
1881             return data_;
1882         }
1883         BOOST_UBLAS_INLINE
data()1884         matrix_closure_type &data () {
1885             return data_;
1886         }
1887 
1888         // Element access
1889 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
1890         BOOST_UBLAS_INLINE
operator ()(size_type i) const1891         const_reference operator () (size_type i) const {
1892             return data_ (ia1_ (i), ia2_ (i));
1893         }
1894         BOOST_UBLAS_INLINE
operator ()(size_type i)1895         reference operator () (size_type i) {
1896             return data_ (ia1_ (i), ia2_ (i));
1897         }
1898 
1899         BOOST_UBLAS_INLINE
operator [](size_type i) const1900         const_reference operator [] (size_type i) const {
1901             return (*this) (i);
1902         }
1903         BOOST_UBLAS_INLINE
operator [](size_type i)1904         reference operator [] (size_type i) {
1905             return (*this) (i);
1906         }
1907 #else
1908         BOOST_UBLAS_INLINE
operator ()(size_type i) const1909         reference operator () (size_type i) const {
1910             return data_ (ia1_ (i), ia2_ (i));
1911         }
1912 
1913         BOOST_UBLAS_INLINE
operator [](size_type i) const1914         reference operator [] (size_type i) const {
1915             return (*this) (i);
1916         }
1917 #endif
1918 
1919         // Assignment
1920         BOOST_UBLAS_INLINE
operator =(const matrix_vector_indirect & mvi)1921         matrix_vector_indirect &operator = (const matrix_vector_indirect &mvi) {
1922             // ISSUE need a temporary, proxy can be overlaping alias
1923             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (mvi));
1924             return *this;
1925         }
1926         BOOST_UBLAS_INLINE
assign_temporary(matrix_vector_indirect & mvi)1927         matrix_vector_indirect &assign_temporary (matrix_vector_indirect &mvi) {
1928             // assign elements, proxied container remains the same
1929             vector_assign<scalar_assign> (*this, mvi);
1930             return *this;
1931         }
1932         template<class AE>
1933         BOOST_UBLAS_INLINE
operator =(const vector_expression<AE> & ae)1934         matrix_vector_indirect &operator = (const vector_expression<AE> &ae) {
1935             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (ae));
1936             return *this;
1937         }
1938         template<class AE>
1939         BOOST_UBLAS_INLINE
assign(const vector_expression<AE> & ae)1940         matrix_vector_indirect &assign (const vector_expression<AE> &ae) {
1941             vector_assign<scalar_assign> (*this, ae);
1942             return *this;
1943         }
1944         template<class AE>
1945         BOOST_UBLAS_INLINE
operator +=(const vector_expression<AE> & ae)1946         matrix_vector_indirect &operator += (const vector_expression<AE> &ae) {
1947             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this + ae));
1948             return *this;
1949         }
1950         template<class AE>
1951         BOOST_UBLAS_INLINE
plus_assign(const vector_expression<AE> & ae)1952         matrix_vector_indirect &plus_assign (const vector_expression<AE> &ae) {
1953             vector_assign<scalar_plus_assign> (*this, ae);
1954             return *this;
1955         }
1956         template<class AE>
1957         BOOST_UBLAS_INLINE
operator -=(const vector_expression<AE> & ae)1958         matrix_vector_indirect &operator -= (const vector_expression<AE> &ae) {
1959             vector_assign<scalar_assign> (*this, typename vector_temporary_traits<M>::type (*this - ae));
1960             return *this;
1961         }
1962         template<class AE>
1963         BOOST_UBLAS_INLINE
minus_assign(const vector_expression<AE> & ae)1964         matrix_vector_indirect &minus_assign (const vector_expression<AE> &ae) {
1965             vector_assign<scalar_minus_assign> (*this, ae);
1966             return *this;
1967         }
1968         template<class AT>
1969         BOOST_UBLAS_INLINE
operator *=(const AT & at)1970         matrix_vector_indirect &operator *= (const AT &at) {
1971             vector_assign_scalar<scalar_multiplies_assign> (*this, at);
1972             return *this;
1973         }
1974         template<class AT>
1975         BOOST_UBLAS_INLINE
operator /=(const AT & at)1976         matrix_vector_indirect &operator /= (const AT &at) {
1977             vector_assign_scalar<scalar_divides_assign> (*this, at);
1978             return *this;
1979         }
1980 
1981         // Closure comparison
1982         BOOST_UBLAS_INLINE
same_closure(const matrix_vector_indirect & mvi) const1983         bool same_closure (const matrix_vector_indirect &mvi) const {
1984             return (*this).data_.same_closure (mvi.data_);
1985         }
1986 
1987         // Comparison
1988         BOOST_UBLAS_INLINE
operator ==(const matrix_vector_indirect & mvi) const1989         bool operator == (const matrix_vector_indirect &mvi) const {
1990             return (*this).data_ == mvi.data_ && ia1_ == mvi.ia1_ && ia2_ == mvi.ia2_;
1991         }
1992 
1993         // Swapping
1994         BOOST_UBLAS_INLINE
swap(matrix_vector_indirect mvi)1995         void swap (matrix_vector_indirect mvi) {
1996             if (this != &mvi) {
1997                 BOOST_UBLAS_CHECK (size () == mvi.size (), bad_size ());
1998                 // Sparse ranges may be nonconformant now.
1999                 // std::swap_ranges (begin (), end (), mvi.begin ());
2000                 vector_swap<scalar_swap> (*this, mvi);
2001             }
2002         }
2003         BOOST_UBLAS_INLINE
swap(matrix_vector_indirect mvi1,matrix_vector_indirect mvi2)2004         friend void swap (matrix_vector_indirect mvi1, matrix_vector_indirect mvi2) {
2005             mvi1.swap (mvi2);
2006         }
2007 
2008         // Iterator types
2009     private:
2010         // Use indirect array as an index - FIXME this fails for packed assignment
2011         typedef typename IA::const_iterator const_subiterator1_type;
2012         typedef typename IA::const_iterator subiterator1_type;
2013         typedef typename IA::const_iterator const_subiterator2_type;
2014         typedef typename IA::const_iterator subiterator2_type;
2015 
2016     public:
2017         class const_iterator;
2018         class iterator;
2019 
2020         // Element lookup
2021         BOOST_UBLAS_INLINE
find(size_type i) const2022         const_iterator find (size_type i) const {
2023             return const_iterator (*this, ia1_.begin () + i, ia2_.begin () + i);
2024         }
2025         BOOST_UBLAS_INLINE
find(size_type i)2026         iterator find (size_type i) {
2027             return iterator (*this, ia1_.begin () + i, ia2_.begin () + i);
2028         }
2029 
2030         // Iterators simply are indices.
2031 
2032         class const_iterator:
2033             public container_const_reference<matrix_vector_indirect>,
2034             public iterator_base_traits<typename M::const_iterator1::iterator_category>::template
2035                         iterator_base<const_iterator, value_type>::type {
2036         public:
2037             // FIXME Iterator can never be different code was:
2038             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
2039             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
2040 
2041             typedef typename matrix_vector_indirect::value_type value_type;
2042             typedef typename matrix_vector_indirect::difference_type difference_type;
2043             typedef typename matrix_vector_indirect::const_reference reference;
2044             typedef const typename matrix_vector_indirect::value_type *pointer;
2045 
2046             // Construction and destruction
2047             BOOST_UBLAS_INLINE
const_iterator()2048             const_iterator ():
2049                 container_const_reference<self_type> (), it1_ (), it2_ () {}
2050             BOOST_UBLAS_INLINE
const_iterator(const self_type & mvi,const const_subiterator1_type & it1,const const_subiterator2_type & it2)2051             const_iterator (const self_type &mvi, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
2052                 container_const_reference<self_type> (mvi), it1_ (it1), it2_ (it2) {}
2053             BOOST_UBLAS_INLINE
const_iterator(const iterator & it)2054             const_iterator (const iterator &it):
2055                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
2056 
2057             // Arithmetic
2058             BOOST_UBLAS_INLINE
operator ++()2059             const_iterator &operator ++ () {
2060                 ++ it1_;
2061                 ++ it2_;
2062                 return *this;
2063             }
2064             BOOST_UBLAS_INLINE
operator --()2065             const_iterator &operator -- () {
2066                 -- it1_;
2067                 -- it2_;
2068                 return *this;
2069             }
2070             BOOST_UBLAS_INLINE
operator +=(difference_type n)2071             const_iterator &operator += (difference_type n) {
2072                 it1_ += n;
2073                 it2_ += n;
2074                 return *this;
2075             }
2076             BOOST_UBLAS_INLINE
operator -=(difference_type n)2077             const_iterator &operator -= (difference_type n) {
2078                 it1_ -= n;
2079                 it2_ -= n;
2080                 return *this;
2081             }
2082             BOOST_UBLAS_INLINE
operator -(const const_iterator & it) const2083             difference_type operator - (const const_iterator &it) const {
2084                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2085                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
2086             }
2087 
2088             // Dereference
2089             BOOST_UBLAS_INLINE
operator *() const2090             const_reference operator * () const {
2091                 // FIXME replace find with at_element
2092                 return (*this) ().data_ (*it1_, *it2_);
2093             }
2094 
2095             // Index
2096             BOOST_UBLAS_INLINE
index() const2097             size_type  index () const {
2098                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
2099             }
2100 
2101             // Assignment
2102             BOOST_UBLAS_INLINE
operator =(const const_iterator & it)2103             const_iterator &operator = (const const_iterator &it) {
2104                 container_const_reference<self_type>::assign (&it ());
2105                 it1_ = it.it1_;
2106                 it2_ = it.it2_;
2107                 return *this;
2108             }
2109 
2110             // Comparison
2111             BOOST_UBLAS_INLINE
operator ==(const const_iterator & it) const2112             bool operator == (const const_iterator &it) const {
2113                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2114                 return it1_ == it.it1_ && it2_ == it.it2_;
2115             }
2116             BOOST_UBLAS_INLINE
operator <(const const_iterator & it) const2117             bool operator < (const const_iterator &it) const {
2118                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2119                 return it1_ < it.it1_ && it2_ < it.it2_;
2120             }
2121 
2122         private:
2123             const_subiterator1_type it1_;
2124             const_subiterator2_type it2_;
2125         };
2126 
2127         BOOST_UBLAS_INLINE
begin() const2128         const_iterator begin () const {
2129             return find (0);
2130         }
2131         BOOST_UBLAS_INLINE
end() const2132         const_iterator end () const {
2133             return find (size ());
2134         }
2135 
2136         class iterator:
2137             public container_reference<matrix_vector_indirect>,
2138             public iterator_base_traits<typename M::iterator1::iterator_category>::template
2139                         iterator_base<iterator, value_type>::type {
2140         public:
2141             // FIXME Iterator can never be different code was:
2142             // typename iterator_restrict_traits<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::iterator_category>
2143             BOOST_STATIC_ASSERT ((boost::is_same<typename M::const_iterator1::iterator_category, typename M::const_iterator2::iterator_category>::value ));
2144 
2145             typedef typename matrix_vector_indirect::value_type value_type;
2146             typedef typename matrix_vector_indirect::difference_type difference_type;
2147             typedef typename matrix_vector_indirect::reference reference;
2148             typedef typename matrix_vector_indirect::value_type *pointer;
2149 
2150             // Construction and destruction
2151             BOOST_UBLAS_INLINE
iterator()2152             iterator ():
2153                 container_reference<self_type> (), it1_ (), it2_ () {}
2154             BOOST_UBLAS_INLINE
iterator(self_type & mvi,const subiterator1_type & it1,const subiterator2_type & it2)2155             iterator (self_type &mvi, const subiterator1_type &it1, const subiterator2_type &it2):
2156                 container_reference<self_type> (mvi), it1_ (it1), it2_ (it2) {}
2157 
2158             // Arithmetic
2159             BOOST_UBLAS_INLINE
operator ++()2160             iterator &operator ++ () {
2161                 ++ it1_;
2162                 ++ it2_;
2163                 return *this;
2164             }
2165             BOOST_UBLAS_INLINE
operator --()2166             iterator &operator -- () {
2167                 -- it1_;
2168                 -- it2_;
2169                 return *this;
2170             }
2171             BOOST_UBLAS_INLINE
operator +=(difference_type n)2172             iterator &operator += (difference_type n) {
2173                 it1_ += n;
2174                 it2_ += n;
2175                 return *this;
2176             }
2177             BOOST_UBLAS_INLINE
operator -=(difference_type n)2178             iterator &operator -= (difference_type n) {
2179                 it1_ -= n;
2180                 it2_ -= n;
2181                 return *this;
2182             }
2183             BOOST_UBLAS_INLINE
operator -(const iterator & it) const2184             difference_type operator - (const iterator &it) const {
2185                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2186                 return BOOST_UBLAS_SAME (it1_ - it.it1_, it2_ - it.it2_);
2187             }
2188 
2189             // Dereference
2190             BOOST_UBLAS_INLINE
operator *() const2191             reference operator * () const {
2192                 // FIXME replace find with at_element
2193                 return (*this) ().data_ (*it1_, *it2_);
2194             }
2195 
2196             // Index
2197             BOOST_UBLAS_INLINE
index() const2198             size_type index () const {
2199                 return BOOST_UBLAS_SAME (it1_.index (), it2_.index ());
2200             }
2201 
2202             // Assignment
2203             BOOST_UBLAS_INLINE
operator =(const iterator & it)2204             iterator &operator = (const iterator &it) {
2205                 container_reference<self_type>::assign (&it ());
2206                 it1_ = it.it1_;
2207                 it2_ = it.it2_;
2208                 return *this;
2209             }
2210 
2211             // Comparison
2212             BOOST_UBLAS_INLINE
operator ==(const iterator & it) const2213             bool operator == (const iterator &it) const {
2214                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2215                 return it1_ == it.it1_ && it2_ == it.it2_;
2216             }
2217             BOOST_UBLAS_INLINE
operator <(const iterator & it) const2218             bool operator < (const iterator &it) const {
2219                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2220                 return it1_ < it.it1_ && it2_ < it.it2_;
2221             }
2222 
2223         private:
2224             subiterator1_type it1_;
2225             subiterator2_type it2_;
2226 
2227             friend class const_iterator;
2228         };
2229 
2230         BOOST_UBLAS_INLINE
begin()2231         iterator begin () {
2232             return find (0);
2233         }
2234         BOOST_UBLAS_INLINE
end()2235         iterator end () {
2236             return find (size ());
2237         }
2238 
2239         // Reverse iterator
2240         typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
2241         typedef reverse_iterator_base<iterator> reverse_iterator;
2242 
2243         BOOST_UBLAS_INLINE
rbegin() const2244         const_reverse_iterator rbegin () const {
2245             return const_reverse_iterator (end ());
2246         }
2247         BOOST_UBLAS_INLINE
rend() const2248         const_reverse_iterator rend () const {
2249             return const_reverse_iterator (begin ());
2250         }
2251         BOOST_UBLAS_INLINE
rbegin()2252         reverse_iterator rbegin () {
2253             return reverse_iterator (end ());
2254         }
2255         BOOST_UBLAS_INLINE
rend()2256         reverse_iterator rend () {
2257             return reverse_iterator (begin ());
2258         }
2259 
2260     private:
2261         matrix_closure_type data_;
2262         indirect_array_type ia1_;
2263         indirect_array_type ia2_;
2264     };
2265 
2266     // Specialize temporary
2267     template <class M, class IA>
2268     struct vector_temporary_traits< matrix_vector_indirect<M,IA> >
2269     : vector_temporary_traits< M > {} ;
2270 
2271     // Matrix based range class
2272     template<class M>
2273     class matrix_range:
2274         public matrix_expression<matrix_range<M> > {
2275 
2276         typedef matrix_range<M> self_type;
2277     public:
2278 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
2279         using matrix_expression<self_type>::operator ();
2280 #endif
2281         typedef M matrix_type;
2282         typedef typename M::size_type size_type;
2283         typedef typename M::difference_type difference_type;
2284         typedef typename M::value_type value_type;
2285         typedef typename M::const_reference const_reference;
2286         typedef typename boost::mpl::if_<boost::is_const<M>,
2287                                           typename M::const_reference,
2288                                           typename M::reference>::type reference;
2289         typedef typename boost::mpl::if_<boost::is_const<M>,
2290                                           typename M::const_closure_type,
2291                                           typename M::closure_type>::type matrix_closure_type;
2292         typedef basic_range<size_type, difference_type> range_type;
2293         typedef const self_type const_closure_type;
2294         typedef self_type closure_type;
2295         typedef typename storage_restrict_traits<typename M::storage_category,
2296                                                  dense_proxy_tag>::storage_category storage_category;
2297         typedef typename M::orientation_category orientation_category;
2298 
2299         // Construction and destruction
2300         BOOST_UBLAS_INLINE
matrix_range(matrix_type & data,const range_type & r1,const range_type & r2)2301         matrix_range (matrix_type &data, const range_type &r1, const range_type &r2):
2302             data_ (data), r1_ (r1.preprocess (data.size1 ())), r2_ (r2.preprocess (data.size2 ())) {
2303             // Early checking of preconditions here.
2304             // BOOST_UBLAS_CHECK (r1_.start () <= data_.size1 () &&
2305             //                    r1_.start () + r1_.size () <= data_.size1 (), bad_index ());
2306             // BOOST_UBLAS_CHECK (r2_.start () <= data_.size2 () &&
2307             //                    r2_.start () + r2_.size () <= data_.size2 (), bad_index ());
2308         }
2309         BOOST_UBLAS_INLINE
matrix_range(const matrix_closure_type & data,const range_type & r1,const range_type & r2,int)2310         matrix_range (const matrix_closure_type &data, const range_type &r1, const range_type &r2, int):
2311             data_ (data), r1_ (r1.preprocess (data.size1 ())), r2_ (r2.preprocess (data.size2 ())) {
2312             // Early checking of preconditions here.
2313             // BOOST_UBLAS_CHECK (r1_.start () <= data_.size1 () &&
2314             //                    r1_.start () + r1_.size () <= data_.size1 (), bad_index ());
2315             // BOOST_UBLAS_CHECK (r2_.start () <= data_.size2 () &&
2316             //                    r2_.start () + r2_.size () <= data_.size2 (), bad_index ());
2317         }
2318 
2319         // Accessors
2320         BOOST_UBLAS_INLINE
start1() const2321         size_type start1 () const {
2322             return r1_.start ();
2323         }
2324         BOOST_UBLAS_INLINE
size1() const2325         size_type size1 () const {
2326             return r1_.size ();
2327         }
2328         BOOST_UBLAS_INLINE
start2() const2329         size_type start2() const {
2330             return r2_.start ();
2331         }
2332         BOOST_UBLAS_INLINE
size2() const2333         size_type size2 () const {
2334             return r2_.size ();
2335         }
2336 
2337         // Storage accessors
2338         BOOST_UBLAS_INLINE
data() const2339         const matrix_closure_type &data () const {
2340             return data_;
2341         }
2342         BOOST_UBLAS_INLINE
data()2343         matrix_closure_type &data () {
2344             return data_;
2345         }
2346 
2347         // Element access
2348 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
2349         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const2350         const_reference operator () (size_type i, size_type j) const {
2351             return data_ (r1_ (i), r2_ (j));
2352         }
2353         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j)2354         reference operator () (size_type i, size_type j) {
2355             return data_ (r1_ (i), r2_ (j));
2356         }
2357 #else
2358         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const2359         reference operator () (size_type i, size_type j) const {
2360             return data_ (r1_ (i), r2_ (j));
2361         }
2362 #endif
2363 
2364         // ISSUE can this be done in free project function?
2365         // Although a const function can create a non-const proxy to a non-const object
2366         // Critical is that matrix_type and data_ (vector_closure_type) are const correct
2367         BOOST_UBLAS_INLINE
project(const range_type & r1,const range_type & r2) const2368         matrix_range<matrix_type> project (const range_type &r1, const range_type &r2) const {
2369             return matrix_range<matrix_type>  (data_, r1_.compose (r1.preprocess (data_.size1 ())), r2_.compose (r2.preprocess (data_.size2 ())), 0);
2370         }
2371 
2372         // Assignment
2373         BOOST_UBLAS_INLINE
operator =(const matrix_range & mr)2374         matrix_range &operator = (const matrix_range &mr) {
2375             matrix_assign<scalar_assign> (*this, mr);
2376             return *this;
2377         }
2378         BOOST_UBLAS_INLINE
assign_temporary(matrix_range & mr)2379         matrix_range &assign_temporary (matrix_range &mr) {
2380             return *this = mr;
2381         }
2382         template<class AE>
2383         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)2384         matrix_range &operator = (const matrix_expression<AE> &ae) {
2385             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (ae));
2386             return *this;
2387         }
2388         template<class AE>
2389         BOOST_UBLAS_INLINE
assign(const matrix_expression<AE> & ae)2390         matrix_range &assign (const matrix_expression<AE> &ae) {
2391             matrix_assign<scalar_assign> (*this, ae);
2392             return *this;
2393         }
2394         template<class AE>
2395         BOOST_UBLAS_INLINE
operator +=(const matrix_expression<AE> & ae)2396         matrix_range& operator += (const matrix_expression<AE> &ae) {
2397             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this + ae));
2398             return *this;
2399         }
2400         template<class AE>
2401         BOOST_UBLAS_INLINE
plus_assign(const matrix_expression<AE> & ae)2402         matrix_range &plus_assign (const matrix_expression<AE> &ae) {
2403             matrix_assign<scalar_plus_assign> (*this, ae);
2404             return *this;
2405         }
2406         template<class AE>
2407         BOOST_UBLAS_INLINE
operator -=(const matrix_expression<AE> & ae)2408         matrix_range& operator -= (const matrix_expression<AE> &ae) {
2409             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this - ae));
2410             return *this;
2411         }
2412         template<class AE>
2413         BOOST_UBLAS_INLINE
minus_assign(const matrix_expression<AE> & ae)2414         matrix_range &minus_assign (const matrix_expression<AE> &ae) {
2415             matrix_assign<scalar_minus_assign> (*this, ae);
2416             return *this;
2417         }
2418         template<class AT>
2419         BOOST_UBLAS_INLINE
operator *=(const AT & at)2420         matrix_range& operator *= (const AT &at) {
2421             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
2422             return *this;
2423         }
2424         template<class AT>
2425         BOOST_UBLAS_INLINE
operator /=(const AT & at)2426         matrix_range& operator /= (const AT &at) {
2427             matrix_assign_scalar<scalar_divides_assign> (*this, at);
2428             return *this;
2429         }
2430 
2431         // Closure comparison
2432         BOOST_UBLAS_INLINE
same_closure(const matrix_range & mr) const2433         bool same_closure (const matrix_range &mr) const {
2434             return (*this).data_.same_closure (mr.data_);
2435         }
2436 
2437         // Comparison
2438         BOOST_UBLAS_INLINE
operator ==(const matrix_range & mr) const2439         bool operator == (const matrix_range &mr) const {
2440             return (*this).data_ == (mr.data_) && r1_ == mr.r1_ && r2_ == mr.r2_;
2441         }
2442 
2443         // Swapping
2444         BOOST_UBLAS_INLINE
swap(matrix_range mr)2445         void swap (matrix_range mr) {
2446             if (this != &mr) {
2447                 BOOST_UBLAS_CHECK (size1 () == mr.size1 (), bad_size ());
2448                 BOOST_UBLAS_CHECK (size2 () == mr.size2 (), bad_size ());
2449                 matrix_swap<scalar_swap> (*this, mr);
2450             }
2451         }
2452         BOOST_UBLAS_INLINE
swap(matrix_range mr1,matrix_range mr2)2453         friend void swap (matrix_range mr1, matrix_range mr2) {
2454             mr1.swap (mr2);
2455         }
2456 
2457         // Iterator types
2458     private:
2459         typedef typename M::const_iterator1 const_subiterator1_type;
2460         typedef typename boost::mpl::if_<boost::is_const<M>,
2461                                           typename M::const_iterator1,
2462                                           typename M::iterator1>::type subiterator1_type;
2463         typedef typename M::const_iterator2 const_subiterator2_type;
2464         typedef typename boost::mpl::if_<boost::is_const<M>,
2465                                           typename M::const_iterator2,
2466                                           typename M::iterator2>::type subiterator2_type;
2467 
2468     public:
2469 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
2470         typedef indexed_iterator1<matrix_range<matrix_type>,
2471                                   typename subiterator1_type::iterator_category> iterator1;
2472         typedef indexed_iterator2<matrix_range<matrix_type>,
2473                                   typename subiterator2_type::iterator_category> iterator2;
2474         typedef indexed_const_iterator1<matrix_range<matrix_type>,
2475                                         typename const_subiterator1_type::iterator_category> const_iterator1;
2476         typedef indexed_const_iterator2<matrix_range<matrix_type>,
2477                                         typename const_subiterator2_type::iterator_category> const_iterator2;
2478 #else
2479         class const_iterator1;
2480         class iterator1;
2481         class const_iterator2;
2482         class iterator2;
2483 #endif
2484         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
2485         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
2486         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
2487         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
2488 
2489         // Element lookup
2490         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j) const2491         const_iterator1 find1 (int rank, size_type i, size_type j) const {
2492             const_subiterator1_type it1 (data_.find1 (rank, start1 () + i, start2 () + j));
2493 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
2494             return const_iterator1 (*this, it1.index1 (), it1.index2 ());
2495 #else
2496             return const_iterator1 (*this, it1);
2497 #endif
2498         }
2499         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j)2500         iterator1 find1 (int rank, size_type i, size_type j) {
2501             subiterator1_type it1 (data_.find1 (rank, start1 () + i, start2 () + j));
2502 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
2503             return iterator1 (*this, it1.index1 (), it1.index2 ());
2504 #else
2505             return iterator1 (*this, it1);
2506 #endif
2507         }
2508         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j) const2509         const_iterator2 find2 (int rank, size_type i, size_type j) const {
2510             const_subiterator2_type it2 (data_.find2 (rank, start1 () + i, start2 () + j));
2511 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
2512             return const_iterator2 (*this, it2.index1 (), it2.index2 ());
2513 #else
2514             return const_iterator2 (*this, it2);
2515 #endif
2516         }
2517         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j)2518         iterator2 find2 (int rank, size_type i, size_type j) {
2519             subiterator2_type it2 (data_.find2 (rank, start1 () + i, start2 () + j));
2520 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
2521             return iterator2 (*this, it2.index1 (), it2.index2 ());
2522 #else
2523             return iterator2 (*this, it2);
2524 #endif
2525         }
2526 
2527 
2528 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2529         class const_iterator1:
2530             public container_const_reference<matrix_range>,
2531             public iterator_base_traits<typename const_subiterator1_type::iterator_category>::template
2532                         iterator_base<const_iterator1, value_type>::type {
2533         public:
2534             typedef typename const_subiterator1_type::value_type value_type;
2535             typedef typename const_subiterator1_type::difference_type difference_type;
2536             typedef typename const_subiterator1_type::reference reference;
2537             typedef typename const_subiterator1_type::pointer pointer;
2538             typedef const_iterator2 dual_iterator_type;
2539             typedef const_reverse_iterator2 dual_reverse_iterator_type;
2540 
2541             // Construction and destruction
2542             BOOST_UBLAS_INLINE
const_iterator1()2543             const_iterator1 ():
2544                 container_const_reference<self_type> (), it_ () {}
2545             BOOST_UBLAS_INLINE
const_iterator1(const self_type & mr,const const_subiterator1_type & it)2546             const_iterator1 (const self_type &mr, const const_subiterator1_type &it):
2547                 container_const_reference<self_type> (mr), it_ (it) {}
2548             BOOST_UBLAS_INLINE
const_iterator1(const iterator1 & it)2549             const_iterator1 (const iterator1 &it):
2550                 container_const_reference<self_type> (it ()), it_ (it.it_) {}
2551 
2552             // Arithmetic
2553             BOOST_UBLAS_INLINE
operator ++()2554             const_iterator1 &operator ++ () {
2555                 ++ it_;
2556                 return *this;
2557             }
2558             BOOST_UBLAS_INLINE
operator --()2559             const_iterator1 &operator -- () {
2560                 -- it_;
2561                 return *this;
2562             }
2563             BOOST_UBLAS_INLINE
operator +=(difference_type n)2564             const_iterator1 &operator += (difference_type n) {
2565                 it_ += n;
2566                 return *this;
2567             }
2568             BOOST_UBLAS_INLINE
operator -=(difference_type n)2569             const_iterator1 &operator -= (difference_type n) {
2570                 it_ -= n;
2571                 return *this;
2572             }
2573             BOOST_UBLAS_INLINE
operator -(const const_iterator1 & it) const2574             difference_type operator - (const const_iterator1 &it) const {
2575                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2576                 return it_ - it.it_;
2577             }
2578 
2579             // Dereference
2580             BOOST_UBLAS_INLINE
operator *() const2581             const_reference operator * () const {
2582                 return *it_;
2583             }
2584 
2585 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2586             BOOST_UBLAS_INLINE
2587 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2588             typename self_type::
2589 #endif
begin() const2590             const_iterator2 begin () const {
2591                 const self_type &mr = (*this) ();
2592                 return mr.find2 (1, index1 (), 0);
2593             }
2594             BOOST_UBLAS_INLINE
2595 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2596             typename self_type::
2597 #endif
end() const2598             const_iterator2 end () const {
2599                 const self_type &mr = (*this) ();
2600                 return mr.find2 (1, index1 (), mr.size2 ());
2601             }
2602             BOOST_UBLAS_INLINE
2603 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2604             typename self_type::
2605 #endif
rbegin() const2606             const_reverse_iterator2 rbegin () const {
2607                 return const_reverse_iterator2 (end ());
2608             }
2609             BOOST_UBLAS_INLINE
2610 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2611             typename self_type::
2612 #endif
rend() const2613             const_reverse_iterator2 rend () const {
2614                 return const_reverse_iterator2 (begin ());
2615             }
2616 #endif
2617 
2618             // Indices
2619             BOOST_UBLAS_INLINE
index1() const2620             size_type index1 () const {
2621                 return it_.index1 () - (*this) ().start1 ();
2622             }
2623             BOOST_UBLAS_INLINE
index2() const2624             size_type index2 () const {
2625                 return it_.index2 () - (*this) ().start2 ();
2626             }
2627 
2628             // Assignment
2629             BOOST_UBLAS_INLINE
operator =(const const_iterator1 & it)2630             const_iterator1 &operator = (const const_iterator1 &it) {
2631                 container_const_reference<self_type>::assign (&it ());
2632                 it_ = it.it_;
2633                 return *this;
2634             }
2635 
2636             // Comparison
2637             BOOST_UBLAS_INLINE
operator ==(const const_iterator1 & it) const2638             bool operator == (const const_iterator1 &it) const {
2639                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2640                 return it_ == it.it_;
2641             }
2642             BOOST_UBLAS_INLINE
operator <(const const_iterator1 & it) const2643             bool operator < (const const_iterator1 &it) const {
2644                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2645                 return it_ < it.it_;
2646             }
2647 
2648         private:
2649             const_subiterator1_type it_;
2650         };
2651 #endif
2652 
2653         BOOST_UBLAS_INLINE
begin1() const2654         const_iterator1 begin1 () const {
2655             return find1 (0, 0, 0);
2656         }
2657         BOOST_UBLAS_INLINE
end1() const2658         const_iterator1 end1 () const {
2659             return find1 (0, size1 (), 0);
2660         }
2661 
2662 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2663         class iterator1:
2664             public container_reference<matrix_range>,
2665             public iterator_base_traits<typename subiterator1_type::iterator_category>::template
2666                         iterator_base<iterator1, value_type>::type {
2667         public:
2668             typedef typename subiterator1_type::value_type value_type;
2669             typedef typename subiterator1_type::difference_type difference_type;
2670             typedef typename subiterator1_type::reference reference;
2671             typedef typename subiterator1_type::pointer pointer;
2672             typedef iterator2 dual_iterator_type;
2673             typedef reverse_iterator2 dual_reverse_iterator_type;
2674 
2675             // Construction and destruction
2676             BOOST_UBLAS_INLINE
iterator1()2677             iterator1 ():
2678                 container_reference<self_type> (), it_ () {}
2679             BOOST_UBLAS_INLINE
iterator1(self_type & mr,const subiterator1_type & it)2680             iterator1 (self_type &mr, const subiterator1_type &it):
2681                 container_reference<self_type> (mr), it_ (it) {}
2682 
2683             // Arithmetic
2684             BOOST_UBLAS_INLINE
operator ++()2685             iterator1 &operator ++ () {
2686                 ++ it_;
2687                 return *this;
2688             }
2689             BOOST_UBLAS_INLINE
operator --()2690             iterator1 &operator -- () {
2691                 -- it_;
2692                 return *this;
2693             }
2694             BOOST_UBLAS_INLINE
operator +=(difference_type n)2695             iterator1 &operator += (difference_type n) {
2696                 it_ += n;
2697                 return *this;
2698             }
2699             BOOST_UBLAS_INLINE
operator -=(difference_type n)2700             iterator1 &operator -= (difference_type n) {
2701                 it_ -= n;
2702                 return *this;
2703             }
2704             BOOST_UBLAS_INLINE
operator -(const iterator1 & it) const2705             difference_type operator - (const iterator1 &it) const {
2706                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2707                 return it_ - it.it_;
2708             }
2709 
2710             // Dereference
2711             BOOST_UBLAS_INLINE
operator *() const2712             reference operator * () const {
2713                 return *it_;
2714             }
2715 
2716 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2717             BOOST_UBLAS_INLINE
2718 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2719             typename self_type::
2720 #endif
begin() const2721             iterator2 begin () const {
2722                 self_type &mr = (*this) ();
2723                 return mr.find2 (1, index1 (), 0);
2724             }
2725             BOOST_UBLAS_INLINE
2726 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2727             typename self_type::
2728 #endif
end() const2729             iterator2 end () const {
2730                 self_type &mr = (*this) ();
2731                 return mr.find2 (1, index1 (), mr.size2 ());
2732             }
2733             BOOST_UBLAS_INLINE
2734 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2735             typename self_type::
2736 #endif
rbegin() const2737             reverse_iterator2 rbegin () const {
2738                 return reverse_iterator2 (end ());
2739             }
2740             BOOST_UBLAS_INLINE
2741 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2742             typename self_type::
2743 #endif
rend() const2744             reverse_iterator2 rend () const {
2745                 return reverse_iterator2 (begin ());
2746             }
2747 #endif
2748 
2749             // Indices
2750             BOOST_UBLAS_INLINE
index1() const2751             size_type index1 () const {
2752                 return it_.index1 () - (*this) ().start1 ();
2753             }
2754             BOOST_UBLAS_INLINE
index2() const2755             size_type index2 () const {
2756                 return it_.index2 () - (*this) ().start2 ();
2757             }
2758 
2759             // Assignment
2760             BOOST_UBLAS_INLINE
operator =(const iterator1 & it)2761             iterator1 &operator = (const iterator1 &it) {
2762                 container_reference<self_type>::assign (&it ());
2763                 it_ = it.it_;
2764                 return *this;
2765             }
2766 
2767             // Comparison
2768             BOOST_UBLAS_INLINE
operator ==(const iterator1 & it) const2769             bool operator == (const iterator1 &it) const {
2770                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2771                 return it_ == it.it_;
2772             }
2773             BOOST_UBLAS_INLINE
operator <(const iterator1 & it) const2774             bool operator < (const iterator1 &it) const {
2775                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2776                 return it_ < it.it_;
2777             }
2778 
2779         private:
2780             subiterator1_type it_;
2781 
2782             friend class const_iterator1;
2783         };
2784 #endif
2785 
2786         BOOST_UBLAS_INLINE
begin1()2787         iterator1 begin1 () {
2788             return find1 (0, 0, 0);
2789         }
2790         BOOST_UBLAS_INLINE
end1()2791         iterator1 end1 () {
2792             return find1 (0, size1 (), 0);
2793         }
2794 
2795 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2796         class const_iterator2:
2797             public container_const_reference<matrix_range>,
2798             public iterator_base_traits<typename const_subiterator2_type::iterator_category>::template
2799                         iterator_base<const_iterator2, value_type>::type {
2800         public:
2801             typedef typename const_subiterator2_type::value_type value_type;
2802             typedef typename const_subiterator2_type::difference_type difference_type;
2803             typedef typename const_subiterator2_type::reference reference;
2804             typedef typename const_subiterator2_type::pointer pointer;
2805             typedef const_iterator1 dual_iterator_type;
2806             typedef const_reverse_iterator1 dual_reverse_iterator_type;
2807 
2808             // Construction and destruction
2809             BOOST_UBLAS_INLINE
const_iterator2()2810             const_iterator2 ():
2811                 container_const_reference<self_type> (), it_ () {}
2812             BOOST_UBLAS_INLINE
const_iterator2(const self_type & mr,const const_subiterator2_type & it)2813             const_iterator2 (const self_type &mr, const const_subiterator2_type &it):
2814                 container_const_reference<self_type> (mr), it_ (it) {}
2815             BOOST_UBLAS_INLINE
const_iterator2(const iterator2 & it)2816             const_iterator2 (const iterator2 &it):
2817                 container_const_reference<self_type> (it ()), it_ (it.it_) {}
2818 
2819             // Arithmetic
2820             BOOST_UBLAS_INLINE
operator ++()2821             const_iterator2 &operator ++ () {
2822                 ++ it_;
2823                 return *this;
2824             }
2825             BOOST_UBLAS_INLINE
operator --()2826             const_iterator2 &operator -- () {
2827                 -- it_;
2828                 return *this;
2829             }
2830             BOOST_UBLAS_INLINE
operator +=(difference_type n)2831             const_iterator2 &operator += (difference_type n) {
2832                 it_ += n;
2833                 return *this;
2834             }
2835             BOOST_UBLAS_INLINE
operator -=(difference_type n)2836             const_iterator2 &operator -= (difference_type n) {
2837                 it_ -= n;
2838                 return *this;
2839             }
2840             BOOST_UBLAS_INLINE
operator -(const const_iterator2 & it) const2841             difference_type operator - (const const_iterator2 &it) const {
2842                 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
2843                 return it_ - it.it_;
2844             }
2845 
2846             // Dereference
2847             BOOST_UBLAS_INLINE
operator *() const2848             const_reference operator * () const {
2849                 return *it_;
2850             }
2851 
2852 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2853             BOOST_UBLAS_INLINE
2854 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2855             typename self_type::
2856 #endif
begin() const2857             const_iterator1 begin () const {
2858                 const self_type &mr = (*this) ();
2859                 return mr.find1 (1, 0, index2 ());
2860             }
2861             BOOST_UBLAS_INLINE
2862 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2863             typename self_type::
2864 #endif
end() const2865             const_iterator1 end () const {
2866                 const self_type &mr = (*this) ();
2867                 return mr.find1 (1, mr.size1 (), index2 ());
2868             }
2869             BOOST_UBLAS_INLINE
2870 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2871             typename self_type::
2872 #endif
rbegin() const2873             const_reverse_iterator1 rbegin () const {
2874                 return const_reverse_iterator1 (end ());
2875             }
2876             BOOST_UBLAS_INLINE
2877 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2878             typename self_type::
2879 #endif
rend() const2880             const_reverse_iterator1 rend () const {
2881                 return const_reverse_iterator1 (begin ());
2882             }
2883 #endif
2884 
2885             // Indices
2886             BOOST_UBLAS_INLINE
index1() const2887             size_type index1 () const {
2888                 return it_.index1 () - (*this) ().start1 ();
2889             }
2890             BOOST_UBLAS_INLINE
index2() const2891             size_type index2 () const {
2892                 return it_.index2 () - (*this) ().start2 ();
2893             }
2894 
2895             // Assignment
2896             BOOST_UBLAS_INLINE
operator =(const const_iterator2 & it)2897             const_iterator2 &operator = (const const_iterator2 &it) {
2898                 container_const_reference<self_type>::assign (&it ());
2899                 it_ = it.it_;
2900                 return *this;
2901             }
2902 
2903             // Comparison
2904             BOOST_UBLAS_INLINE
operator ==(const const_iterator2 & it) const2905             bool operator == (const const_iterator2 &it) const {
2906                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2907                 return it_ == it.it_;
2908             }
2909             BOOST_UBLAS_INLINE
operator <(const const_iterator2 & it) const2910             bool operator < (const const_iterator2 &it) const {
2911                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2912                 return it_ < it.it_;
2913             }
2914 
2915         private:
2916             const_subiterator2_type it_;
2917         };
2918 #endif
2919 
2920         BOOST_UBLAS_INLINE
begin2() const2921         const_iterator2 begin2 () const {
2922             return find2 (0, 0, 0);
2923         }
2924         BOOST_UBLAS_INLINE
end2() const2925         const_iterator2 end2 () const {
2926             return find2 (0, 0, size2 ());
2927         }
2928 
2929 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2930         class iterator2:
2931             public container_reference<matrix_range>,
2932             public iterator_base_traits<typename subiterator2_type::iterator_category>::template
2933                         iterator_base<iterator2, value_type>::type {
2934         public:
2935             typedef typename subiterator2_type::value_type value_type;
2936             typedef typename subiterator2_type::difference_type difference_type;
2937             typedef typename subiterator2_type::reference reference;
2938             typedef typename subiterator2_type::pointer pointer;
2939             typedef iterator1 dual_iterator_type;
2940             typedef reverse_iterator1 dual_reverse_iterator_type;
2941 
2942             // Construction and destruction
2943             BOOST_UBLAS_INLINE
iterator2()2944             iterator2 ():
2945                 container_reference<self_type> (), it_ () {}
2946             BOOST_UBLAS_INLINE
iterator2(self_type & mr,const subiterator2_type & it)2947             iterator2 (self_type &mr, const subiterator2_type &it):
2948                 container_reference<self_type> (mr), it_ (it) {}
2949 
2950             // Arithmetic
2951             BOOST_UBLAS_INLINE
operator ++()2952             iterator2 &operator ++ () {
2953                 ++ it_;
2954                 return *this;
2955             }
2956             BOOST_UBLAS_INLINE
operator --()2957             iterator2 &operator -- () {
2958                 -- it_;
2959                 return *this;
2960             }
2961             BOOST_UBLAS_INLINE
operator +=(difference_type n)2962             iterator2 &operator += (difference_type n) {
2963                 it_ += n;
2964                 return *this;
2965             }
2966             BOOST_UBLAS_INLINE
operator -=(difference_type n)2967             iterator2 &operator -= (difference_type n) {
2968                 it_ -= n;
2969                 return *this;
2970             }
2971             BOOST_UBLAS_INLINE
operator -(const iterator2 & it) const2972             difference_type operator - (const iterator2 &it) const {
2973                BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
2974                 return it_ - it.it_;
2975             }
2976 
2977             // Dereference
2978             BOOST_UBLAS_INLINE
operator *() const2979             reference operator * () const {
2980                 return *it_;
2981             }
2982 
2983 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2984             BOOST_UBLAS_INLINE
2985 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2986             typename self_type::
2987 #endif
begin() const2988             iterator1 begin () const {
2989                 self_type &mr = (*this) ();
2990                 return mr.find1 (1, 0, index2 ());
2991             }
2992             BOOST_UBLAS_INLINE
2993 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2994             typename self_type::
2995 #endif
end() const2996             iterator1 end () const {
2997                 self_type &mr = (*this) ();
2998                 return mr.find1 (1, mr.size1 (), index2 ());
2999             }
3000             BOOST_UBLAS_INLINE
3001 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3002             typename self_type::
3003 #endif
rbegin() const3004             reverse_iterator1 rbegin () const {
3005                 return reverse_iterator1 (end ());
3006             }
3007             BOOST_UBLAS_INLINE
3008 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3009             typename self_type::
3010 #endif
rend() const3011             reverse_iterator1 rend () const {
3012                 return reverse_iterator1 (begin ());
3013             }
3014 #endif
3015 
3016             // Indices
3017             BOOST_UBLAS_INLINE
index1() const3018             size_type index1 () const {
3019                 return it_.index1 () - (*this) ().start1 ();
3020             }
3021             BOOST_UBLAS_INLINE
index2() const3022             size_type index2 () const {
3023                 return it_.index2 () - (*this) ().start2 ();
3024             }
3025 
3026             // Assignment
3027             BOOST_UBLAS_INLINE
operator =(const iterator2 & it)3028             iterator2 &operator = (const iterator2 &it) {
3029                 container_reference<self_type>::assign (&it ());
3030                 it_ = it.it_;
3031                 return *this;
3032             }
3033 
3034             // Comparison
3035             BOOST_UBLAS_INLINE
operator ==(const iterator2 & it) const3036             bool operator == (const iterator2 &it) const {
3037                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3038                 return it_ == it.it_;
3039             }
3040             BOOST_UBLAS_INLINE
operator <(const iterator2 & it) const3041             bool operator < (const iterator2 &it) const {
3042                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3043                 return it_ < it.it_;
3044             }
3045 
3046         private:
3047             subiterator2_type it_;
3048 
3049             friend class const_iterator2;
3050         };
3051 #endif
3052 
3053         BOOST_UBLAS_INLINE
begin2()3054         iterator2 begin2 () {
3055             return find2 (0, 0, 0);
3056         }
3057         BOOST_UBLAS_INLINE
end2()3058         iterator2 end2 () {
3059             return find2 (0, 0, size2 ());
3060         }
3061 
3062         // Reverse iterators
3063 
3064         BOOST_UBLAS_INLINE
rbegin1() const3065         const_reverse_iterator1 rbegin1 () const {
3066             return const_reverse_iterator1 (end1 ());
3067         }
3068         BOOST_UBLAS_INLINE
rend1() const3069         const_reverse_iterator1 rend1 () const {
3070             return const_reverse_iterator1 (begin1 ());
3071         }
3072 
3073         BOOST_UBLAS_INLINE
rbegin1()3074         reverse_iterator1 rbegin1 () {
3075             return reverse_iterator1 (end1 ());
3076         }
3077         BOOST_UBLAS_INLINE
rend1()3078         reverse_iterator1 rend1 () {
3079             return reverse_iterator1 (begin1 ());
3080         }
3081 
3082         BOOST_UBLAS_INLINE
rbegin2() const3083         const_reverse_iterator2 rbegin2 () const {
3084             return const_reverse_iterator2 (end2 ());
3085         }
3086         BOOST_UBLAS_INLINE
rend2() const3087         const_reverse_iterator2 rend2 () const {
3088             return const_reverse_iterator2 (begin2 ());
3089         }
3090 
3091         BOOST_UBLAS_INLINE
rbegin2()3092         reverse_iterator2 rbegin2 () {
3093             return reverse_iterator2 (end2 ());
3094         }
3095         BOOST_UBLAS_INLINE
rend2()3096         reverse_iterator2 rend2 () {
3097             return reverse_iterator2 (begin2 ());
3098         }
3099 
3100     private:
3101         matrix_closure_type data_;
3102         range_type r1_;
3103         range_type r2_;
3104     };
3105 
3106     // Simple Projections
3107     template<class M>
3108     BOOST_UBLAS_INLINE
subrange(M & data,typename M::size_type start1,typename M::size_type stop1,typename M::size_type start2,typename M::size_type stop2)3109     matrix_range<M> subrange (M &data, typename M::size_type start1, typename M::size_type stop1, typename M::size_type start2, typename M::size_type stop2) {
3110         typedef basic_range<typename M::size_type, typename M::difference_type> range_type;
3111         return matrix_range<M> (data, range_type (start1, stop1), range_type (start2, stop2));
3112     }
3113     template<class M>
3114     BOOST_UBLAS_INLINE
subrange(const M & data,typename M::size_type start1,typename M::size_type stop1,typename M::size_type start2,typename M::size_type stop2)3115     matrix_range<const M> subrange (const M &data, typename M::size_type start1, typename M::size_type stop1, typename M::size_type start2, typename M::size_type stop2) {
3116         typedef basic_range<typename M::size_type, typename M::difference_type> range_type;
3117         return matrix_range<const M> (data, range_type (start1, stop1), range_type (start2, stop2));
3118     }
3119 
3120     // Generic Projections
3121     template<class M>
3122     BOOST_UBLAS_INLINE
project(M & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)3123     matrix_range<M> project (M &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
3124         return matrix_range<M> (data, r1, r2);
3125     }
3126     template<class M>
3127     BOOST_UBLAS_INLINE
project(const M & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)3128     const matrix_range<const M> project (const M &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
3129         // ISSUE was: return matrix_range<M> (const_cast<M &> (data), r1, r2);
3130         return matrix_range<const M> (data, r1, r2);
3131     }
3132     template<class M>
3133     BOOST_UBLAS_INLINE
project(matrix_range<M> & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)3134     matrix_range<M> project (matrix_range<M> &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
3135         return data.project (r1, r2);
3136     }
3137     template<class M>
3138     BOOST_UBLAS_INLINE
project(const matrix_range<M> & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)3139     const matrix_range<M> project (const matrix_range<M> &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
3140         return data.project (r1, r2);
3141     }
3142 
3143     // Specialization of temporary_traits
3144     template <class M>
3145     struct matrix_temporary_traits< matrix_range<M> >
3146     : matrix_temporary_traits< M > {} ;
3147 
3148     template <class M>
3149     struct vector_temporary_traits< matrix_range<M> >
3150     : vector_temporary_traits< M > {} ;
3151 
3152     // Matrix based slice class
3153     template<class M>
3154     class matrix_slice:
3155         public matrix_expression<matrix_slice<M> > {
3156 
3157         typedef matrix_slice<M> self_type;
3158     public:
3159 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
3160         using matrix_expression<self_type>::operator ();
3161 #endif
3162         typedef M matrix_type;
3163         typedef typename M::size_type size_type;
3164         typedef typename M::difference_type difference_type;
3165         typedef typename M::value_type value_type;
3166         typedef typename M::const_reference const_reference;
3167         typedef typename boost::mpl::if_<boost::is_const<M>,
3168                                           typename M::const_reference,
3169                                           typename M::reference>::type reference;
3170         typedef typename boost::mpl::if_<boost::is_const<M>,
3171                                           typename M::const_closure_type,
3172                                           typename M::closure_type>::type matrix_closure_type;
3173         typedef basic_range<size_type, difference_type> range_type;
3174         typedef basic_slice<size_type, difference_type> slice_type;
3175         typedef const self_type const_closure_type;
3176         typedef self_type closure_type;
3177         typedef typename storage_restrict_traits<typename M::storage_category,
3178                                                  dense_proxy_tag>::storage_category storage_category;
3179         typedef typename M::orientation_category orientation_category;
3180 
3181         // Construction and destruction
3182         BOOST_UBLAS_INLINE
matrix_slice(matrix_type & data,const slice_type & s1,const slice_type & s2)3183         matrix_slice (matrix_type &data, const slice_type &s1, const slice_type &s2):
3184             data_ (data), s1_ (s1.preprocess (data.size1 ())), s2_ (s2.preprocess (data.size2 ())) {
3185             // Early checking of preconditions here.
3186             // BOOST_UBLAS_CHECK (s1_.start () <= data_.size1 () &&
3187             //                    s1_.start () + s1_.stride () * (s1_.size () - (s1_.size () > 0)) <= data_.size1 (), bad_index ());
3188             // BOOST_UBLAS_CHECK (s2_.start () <= data_.size2 () &&
3189             //                    s2_.start () + s2_.stride () * (s2_.size () - (s2_.size () > 0)) <= data_.size2 (), bad_index ());
3190         }
3191         BOOST_UBLAS_INLINE
matrix_slice(const matrix_closure_type & data,const slice_type & s1,const slice_type & s2,int)3192         matrix_slice (const matrix_closure_type &data, const slice_type &s1, const slice_type &s2, int):
3193             data_ (data), s1_ (s1.preprocess (data.size1 ())), s2_ (s2.preprocess (data.size2 ())) {
3194             // Early checking of preconditions.
3195             // BOOST_UBLAS_CHECK (s1_.start () <= data_.size1 () &&
3196             //                    s1_.start () + s1_.stride () * (s1_.size () - (s1_.size () > 0)) <= data_.size1 (), bad_index ());
3197             // BOOST_UBLAS_CHECK (s2_.start () <= data_.size2 () &&
3198             //                    s2_.start () + s2_.stride () * (s2_.size () - (s2_.size () > 0)) <= data_.size2 (), bad_index ());
3199         }
3200 
3201         // Accessors
3202         BOOST_UBLAS_INLINE
start1() const3203         size_type start1 () const {
3204             return s1_.start ();
3205         }
3206         BOOST_UBLAS_INLINE
start2() const3207         size_type start2 () const {
3208             return s2_.start ();
3209         }
3210         BOOST_UBLAS_INLINE
stride1() const3211         difference_type stride1 () const {
3212             return s1_.stride ();
3213         }
3214         BOOST_UBLAS_INLINE
stride2() const3215         difference_type stride2 () const {
3216             return s2_.stride ();
3217         }
3218         BOOST_UBLAS_INLINE
size1() const3219         size_type size1 () const {
3220             return s1_.size ();
3221         }
3222         BOOST_UBLAS_INLINE
size2() const3223         size_type size2 () const {
3224             return s2_.size ();
3225         }
3226 
3227         // Storage accessors
3228         BOOST_UBLAS_INLINE
data() const3229         const matrix_closure_type &data () const {
3230             return data_;
3231         }
3232         BOOST_UBLAS_INLINE
data()3233         matrix_closure_type &data () {
3234             return data_;
3235         }
3236 
3237         // Element access
3238 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
3239         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const3240         const_reference operator () (size_type i, size_type j) const {
3241             return data_ (s1_ (i), s2_ (j));
3242         }
3243         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j)3244         reference operator () (size_type i, size_type j) {
3245             return data_ (s1_ (i), s2_ (j));
3246         }
3247 #else
3248         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const3249         reference operator () (size_type i, size_type j) const {
3250             return data_ (s1_ (i), s2_ (j));
3251         }
3252 #endif
3253 
3254         // ISSUE can this be done in free project function?
3255         // Although a const function can create a non-const proxy to a non-const object
3256         // Critical is that matrix_type and data_ (vector_closure_type) are const correct
3257         BOOST_UBLAS_INLINE
project(const range_type & r1,const range_type & r2) const3258         matrix_slice<matrix_type> project (const range_type &r1, const range_type &r2) const {
3259             return matrix_slice<matrix_type>  (data_, s1_.compose (r1.preprocess (data_.size1 ())), s2_.compose (r2.preprocess (data_.size2 ())), 0);
3260         }
3261         BOOST_UBLAS_INLINE
project(const slice_type & s1,const slice_type & s2) const3262         matrix_slice<matrix_type> project (const slice_type &s1, const slice_type &s2) const {
3263             return matrix_slice<matrix_type>  (data_, s1_.compose (s1.preprocess (data_.size1 ())), s2_.compose (s2.preprocess (data_.size2 ())), 0);
3264         }
3265 
3266         // Assignment
3267         BOOST_UBLAS_INLINE
operator =(const matrix_slice & ms)3268         matrix_slice &operator = (const matrix_slice &ms) {
3269             matrix_assign<scalar_assign> (*this, ms);
3270             return *this;
3271         }
3272         BOOST_UBLAS_INLINE
assign_temporary(matrix_slice & ms)3273         matrix_slice &assign_temporary (matrix_slice &ms) {
3274             return *this = ms;
3275         }
3276         template<class AE>
3277         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)3278         matrix_slice &operator = (const matrix_expression<AE> &ae) {
3279             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (ae));
3280             return *this;
3281         }
3282         template<class AE>
3283         BOOST_UBLAS_INLINE
assign(const matrix_expression<AE> & ae)3284         matrix_slice &assign (const matrix_expression<AE> &ae) {
3285             matrix_assign<scalar_assign> (*this, ae);
3286             return *this;
3287         }
3288         template<class AE>
3289         BOOST_UBLAS_INLINE
operator +=(const matrix_expression<AE> & ae)3290         matrix_slice& operator += (const matrix_expression<AE> &ae) {
3291             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this + ae));
3292             return *this;
3293         }
3294         template<class AE>
3295         BOOST_UBLAS_INLINE
plus_assign(const matrix_expression<AE> & ae)3296         matrix_slice &plus_assign (const matrix_expression<AE> &ae) {
3297             matrix_assign<scalar_plus_assign> (*this, ae);
3298             return *this;
3299         }
3300         template<class AE>
3301         BOOST_UBLAS_INLINE
operator -=(const matrix_expression<AE> & ae)3302         matrix_slice& operator -= (const matrix_expression<AE> &ae) {
3303             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this - ae));
3304             return *this;
3305         }
3306         template<class AE>
3307         BOOST_UBLAS_INLINE
minus_assign(const matrix_expression<AE> & ae)3308         matrix_slice &minus_assign (const matrix_expression<AE> &ae) {
3309             matrix_assign<scalar_minus_assign> (*this, ae);
3310             return *this;
3311         }
3312         template<class AT>
3313         BOOST_UBLAS_INLINE
operator *=(const AT & at)3314         matrix_slice& operator *= (const AT &at) {
3315             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
3316             return *this;
3317         }
3318         template<class AT>
3319         BOOST_UBLAS_INLINE
operator /=(const AT & at)3320         matrix_slice& operator /= (const AT &at) {
3321             matrix_assign_scalar<scalar_divides_assign> (*this, at);
3322             return *this;
3323         }
3324 
3325         // Closure comparison
3326         BOOST_UBLAS_INLINE
same_closure(const matrix_slice & ms) const3327         bool same_closure (const matrix_slice &ms) const {
3328             return (*this).data_.same_closure (ms.data_);
3329         }
3330 
3331         // Comparison
3332         BOOST_UBLAS_INLINE
operator ==(const matrix_slice & ms) const3333         bool operator == (const matrix_slice &ms) const {
3334             return (*this).data_ == ms.data_ && s1_ == ms.s1_ && s2_ == ms.s2_;
3335         }
3336 
3337         // Swapping
3338         BOOST_UBLAS_INLINE
swap(matrix_slice ms)3339         void swap (matrix_slice ms) {
3340             if (this != &ms) {
3341                 BOOST_UBLAS_CHECK (size1 () == ms.size1 (), bad_size ());
3342                 BOOST_UBLAS_CHECK (size2 () == ms.size2 (), bad_size ());
3343                 matrix_swap<scalar_swap> (*this, ms);
3344             }
3345         }
3346         BOOST_UBLAS_INLINE
swap(matrix_slice ms1,matrix_slice ms2)3347         friend void swap (matrix_slice ms1, matrix_slice ms2) {
3348             ms1.swap (ms2);
3349         }
3350 
3351         // Iterator types
3352     private:
3353         // Use slice as an index - FIXME this fails for packed assignment
3354         typedef typename slice_type::const_iterator const_subiterator1_type;
3355         typedef typename slice_type::const_iterator subiterator1_type;
3356         typedef typename slice_type::const_iterator const_subiterator2_type;
3357         typedef typename slice_type::const_iterator subiterator2_type;
3358 
3359     public:
3360 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
3361         typedef indexed_iterator1<matrix_slice<matrix_type>,
3362                                   typename matrix_type::iterator1::iterator_category> iterator1;
3363         typedef indexed_iterator2<matrix_slice<matrix_type>,
3364                                   typename matrix_type::iterator2::iterator_category> iterator2;
3365         typedef indexed_const_iterator1<matrix_slice<matrix_type>,
3366                                         typename matrix_type::const_iterator1::iterator_category> const_iterator1;
3367         typedef indexed_const_iterator2<matrix_slice<matrix_type>,
3368                                         typename matrix_type::const_iterator2::iterator_category> const_iterator2;
3369 #else
3370         class const_iterator1;
3371         class iterator1;
3372         class const_iterator2;
3373         class iterator2;
3374 #endif
3375         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
3376         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
3377         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
3378         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
3379 
3380         // Element lookup
3381         BOOST_UBLAS_INLINE
find1(int,size_type i,size_type j) const3382         const_iterator1 find1 (int /* rank */, size_type i, size_type j) const {
3383 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
3384             return const_iterator1 (*this, i, j);
3385 #else
3386             return const_iterator1 (*this, s1_.begin () + i, s2_.begin () + j);
3387 #endif
3388         }
3389         BOOST_UBLAS_INLINE
find1(int,size_type i,size_type j)3390         iterator1 find1 (int /* rank */, size_type i, size_type j) {
3391 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
3392             return iterator1 (*this, i, j);
3393 #else
3394             return iterator1 (*this, s1_.begin () + i, s2_.begin () + j);
3395 #endif
3396         }
3397         BOOST_UBLAS_INLINE
find2(int,size_type i,size_type j) const3398         const_iterator2 find2 (int /* rank */, size_type i, size_type j) const {
3399 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
3400             return const_iterator2 (*this, i, j);
3401 #else
3402             return const_iterator2 (*this, s1_.begin () + i, s2_.begin () + j);
3403 #endif
3404         }
3405         BOOST_UBLAS_INLINE
find2(int,size_type i,size_type j)3406         iterator2 find2 (int /* rank */, size_type i, size_type j) {
3407 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
3408             return iterator2 (*this, i, j);
3409 #else
3410             return iterator2 (*this, s1_.begin () + i, s2_.begin () + j);
3411 #endif
3412         }
3413 
3414         // Iterators simply are indices.
3415 
3416 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
3417         class const_iterator1:
3418             public container_const_reference<matrix_slice>,
3419             public iterator_base_traits<typename M::const_iterator1::iterator_category>::template
3420                         iterator_base<const_iterator1, value_type>::type {
3421         public:
3422             typedef typename M::const_iterator1::value_type value_type;
3423             typedef typename M::const_iterator1::difference_type difference_type;
3424             typedef typename M::const_reference reference;    //FIXME due to indexing access
3425             typedef typename M::const_iterator1::pointer pointer;
3426             typedef const_iterator2 dual_iterator_type;
3427             typedef const_reverse_iterator2 dual_reverse_iterator_type;
3428 
3429             // Construction and destruction
3430             BOOST_UBLAS_INLINE
const_iterator1()3431             const_iterator1 ():
3432                 container_const_reference<self_type> (), it1_ (), it2_ () {}
3433             BOOST_UBLAS_INLINE
const_iterator1(const self_type & ms,const const_subiterator1_type & it1,const const_subiterator2_type & it2)3434             const_iterator1 (const self_type &ms, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
3435                 container_const_reference<self_type> (ms), it1_ (it1), it2_ (it2) {}
3436             BOOST_UBLAS_INLINE
const_iterator1(const iterator1 & it)3437             const_iterator1 (const iterator1 &it):
3438                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
3439 
3440             // Arithmetic
3441             BOOST_UBLAS_INLINE
operator ++()3442             const_iterator1 &operator ++ () {
3443                 ++ it1_;
3444                 return *this;
3445             }
3446             BOOST_UBLAS_INLINE
operator --()3447             const_iterator1 &operator -- () {
3448                 -- it1_;
3449                 return *this;
3450             }
3451             BOOST_UBLAS_INLINE
operator +=(difference_type n)3452             const_iterator1 &operator += (difference_type n) {
3453                 it1_ += n;
3454                 return *this;
3455             }
3456             BOOST_UBLAS_INLINE
operator -=(difference_type n)3457             const_iterator1 &operator -= (difference_type n) {
3458                 it1_ -= n;
3459                 return *this;
3460             }
3461             BOOST_UBLAS_INLINE
operator -(const const_iterator1 & it) const3462             difference_type operator - (const const_iterator1 &it) const {
3463                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3464                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3465                 return it1_ - it.it1_;
3466             }
3467 
3468             // Dereference
3469             BOOST_UBLAS_INLINE
operator *() const3470             const_reference operator * () const {
3471                 // FIXME replace find with at_element
3472                 return (*this) ().data_ (*it1_, *it2_);
3473             }
3474 
3475 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
3476             BOOST_UBLAS_INLINE
3477 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3478             typename self_type::
3479 #endif
begin() const3480             const_iterator2 begin () const {
3481                 return const_iterator2 ((*this) (), it1_, it2_ ().begin ());
3482             }
3483             BOOST_UBLAS_INLINE
3484 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3485             typename self_type::
3486 #endif
end() const3487             const_iterator2 end () const {
3488                 return const_iterator2 ((*this) (), it1_, it2_ ().end ());
3489             }
3490             BOOST_UBLAS_INLINE
3491 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3492             typename self_type::
3493 #endif
rbegin() const3494             const_reverse_iterator2 rbegin () const {
3495                 return const_reverse_iterator2 (end ());
3496             }
3497             BOOST_UBLAS_INLINE
3498 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3499             typename self_type::
3500 #endif
rend() const3501             const_reverse_iterator2 rend () const {
3502                 return const_reverse_iterator2 (begin ());
3503             }
3504 #endif
3505 
3506             // Indices
3507             BOOST_UBLAS_INLINE
index1() const3508             size_type index1 () const {
3509                 return it1_.index ();
3510             }
3511             BOOST_UBLAS_INLINE
index2() const3512             size_type index2 () const {
3513                 return it2_.index ();
3514             }
3515 
3516             // Assignment
3517             BOOST_UBLAS_INLINE
operator =(const const_iterator1 & it)3518             const_iterator1 &operator = (const const_iterator1 &it) {
3519                 container_const_reference<self_type>::assign (&it ());
3520                 it1_ = it.it1_;
3521                 it2_ = it.it2_;
3522                 return *this;
3523             }
3524 
3525             // Comparison
3526             BOOST_UBLAS_INLINE
operator ==(const const_iterator1 & it) const3527             bool operator == (const const_iterator1 &it) const {
3528                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3529                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3530                 return it1_ == it.it1_;
3531             }
3532             BOOST_UBLAS_INLINE
operator <(const const_iterator1 & it) const3533             bool operator < (const const_iterator1 &it) const {
3534                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3535                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3536                 return it1_ < it.it1_;
3537             }
3538 
3539         private:
3540             const_subiterator1_type it1_;
3541             const_subiterator2_type it2_;
3542         };
3543 #endif
3544 
3545         BOOST_UBLAS_INLINE
begin1() const3546         const_iterator1 begin1 () const {
3547             return find1 (0, 0, 0);
3548         }
3549         BOOST_UBLAS_INLINE
end1() const3550         const_iterator1 end1 () const {
3551             return find1 (0, size1 (), 0);
3552         }
3553 
3554 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
3555         class iterator1:
3556             public container_reference<matrix_slice>,
3557             public iterator_base_traits<typename M::iterator1::iterator_category>::template
3558                         iterator_base<iterator1, value_type>::type {
3559         public:
3560             typedef typename M::iterator1::value_type value_type;
3561             typedef typename M::iterator1::difference_type difference_type;
3562             typedef typename M::reference reference;    //FIXME due to indexing access
3563             typedef typename M::iterator1::pointer pointer;
3564             typedef iterator2 dual_iterator_type;
3565             typedef reverse_iterator2 dual_reverse_iterator_type;
3566 
3567             // Construction and destruction
3568             BOOST_UBLAS_INLINE
iterator1()3569             iterator1 ():
3570                 container_reference<self_type> (), it1_ (), it2_ () {}
3571             BOOST_UBLAS_INLINE
iterator1(self_type & ms,const subiterator1_type & it1,const subiterator2_type & it2)3572             iterator1 (self_type &ms, const subiterator1_type &it1, const subiterator2_type &it2):
3573                 container_reference<self_type> (ms), it1_ (it1), it2_ (it2) {}
3574 
3575             // Arithmetic
3576             BOOST_UBLAS_INLINE
operator ++()3577             iterator1 &operator ++ () {
3578                 ++ it1_;
3579                 return *this;
3580             }
3581             BOOST_UBLAS_INLINE
operator --()3582             iterator1 &operator -- () {
3583                 -- it1_;
3584                 return *this;
3585             }
3586             BOOST_UBLAS_INLINE
operator +=(difference_type n)3587             iterator1 &operator += (difference_type n) {
3588                 it1_ += n;
3589                 return *this;
3590             }
3591             BOOST_UBLAS_INLINE
operator -=(difference_type n)3592             iterator1 &operator -= (difference_type n) {
3593                 it1_ -= n;
3594                 return *this;
3595             }
3596             BOOST_UBLAS_INLINE
operator -(const iterator1 & it) const3597             difference_type operator - (const iterator1 &it) const {
3598                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3599                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3600                 return it1_ - it.it1_;
3601             }
3602 
3603             // Dereference
3604             BOOST_UBLAS_INLINE
operator *() const3605             reference operator * () const {
3606                 // FIXME replace find with at_element
3607                 return (*this) ().data_ (*it1_, *it2_);
3608             }
3609 
3610 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
3611             BOOST_UBLAS_INLINE
3612 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3613             typename self_type::
3614 #endif
begin() const3615             iterator2 begin () const {
3616                 return iterator2 ((*this) (), it1_, it2_ ().begin ());
3617             }
3618             BOOST_UBLAS_INLINE
3619 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3620             typename self_type::
3621 #endif
end() const3622             iterator2 end () const {
3623                 return iterator2 ((*this) (), it1_, it2_ ().end ());
3624             }
3625             BOOST_UBLAS_INLINE
3626 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3627             typename self_type::
3628 #endif
rbegin() const3629             reverse_iterator2 rbegin () const {
3630                 return reverse_iterator2 (end ());
3631             }
3632             BOOST_UBLAS_INLINE
3633 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3634             typename self_type::
3635 #endif
rend() const3636             reverse_iterator2 rend () const {
3637                 return reverse_iterator2 (begin ());
3638             }
3639 #endif
3640 
3641             // Indices
3642             BOOST_UBLAS_INLINE
index1() const3643             size_type index1 () const {
3644                 return it1_.index ();
3645             }
3646             BOOST_UBLAS_INLINE
index2() const3647             size_type index2 () const {
3648                 return it2_.index ();
3649             }
3650 
3651             // Assignment
3652             BOOST_UBLAS_INLINE
operator =(const iterator1 & it)3653             iterator1 &operator = (const iterator1 &it) {
3654                 container_reference<self_type>::assign (&it ());
3655                 it1_ = it.it1_;
3656                 it2_ = it.it2_;
3657                 return *this;
3658             }
3659 
3660             // Comparison
3661             BOOST_UBLAS_INLINE
operator ==(const iterator1 & it) const3662             bool operator == (const iterator1 &it) const {
3663                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3664                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3665                 return it1_ == it.it1_;
3666             }
3667             BOOST_UBLAS_INLINE
operator <(const iterator1 & it) const3668             bool operator < (const iterator1 &it) const {
3669                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3670                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
3671                 return it1_ < it.it1_;
3672             }
3673 
3674         private:
3675             subiterator1_type it1_;
3676             subiterator2_type it2_;
3677 
3678             friend class const_iterator1;
3679         };
3680 #endif
3681 
3682         BOOST_UBLAS_INLINE
begin1()3683         iterator1 begin1 () {
3684             return find1 (0, 0, 0);
3685         }
3686         BOOST_UBLAS_INLINE
end1()3687         iterator1 end1 () {
3688             return find1 (0, size1 (), 0);
3689         }
3690 
3691 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
3692         class const_iterator2:
3693             public container_const_reference<matrix_slice>,
3694             public iterator_base_traits<typename M::const_iterator2::iterator_category>::template
3695                         iterator_base<const_iterator2, value_type>::type {
3696         public:
3697             typedef typename M::const_iterator2::value_type value_type;
3698             typedef typename M::const_iterator2::difference_type difference_type;
3699             typedef typename M::const_reference reference;    //FIXME due to indexing access
3700             typedef typename M::const_iterator2::pointer pointer;
3701             typedef const_iterator1 dual_iterator_type;
3702             typedef const_reverse_iterator1 dual_reverse_iterator_type;
3703 
3704             // Construction and destruction
3705             BOOST_UBLAS_INLINE
const_iterator2()3706             const_iterator2 ():
3707                 container_const_reference<self_type> (), it1_ (), it2_ () {}
3708             BOOST_UBLAS_INLINE
const_iterator2(const self_type & ms,const const_subiterator1_type & it1,const const_subiterator2_type & it2)3709             const_iterator2 (const self_type &ms, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
3710                 container_const_reference<self_type> (ms), it1_ (it1), it2_ (it2) {}
3711             BOOST_UBLAS_INLINE
const_iterator2(const iterator2 & it)3712             const_iterator2 (const iterator2 &it):
3713                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
3714 
3715             // Arithmetic
3716             BOOST_UBLAS_INLINE
operator ++()3717             const_iterator2 &operator ++ () {
3718                 ++ it2_;
3719                 return *this;
3720             }
3721             BOOST_UBLAS_INLINE
operator --()3722             const_iterator2 &operator -- () {
3723                 -- it2_;
3724                 return *this;
3725             }
3726             BOOST_UBLAS_INLINE
operator +=(difference_type n)3727             const_iterator2 &operator += (difference_type n) {
3728                 it2_ += n;
3729                 return *this;
3730             }
3731             BOOST_UBLAS_INLINE
operator -=(difference_type n)3732             const_iterator2 &operator -= (difference_type n) {
3733                 it2_ -= n;
3734                 return *this;
3735             }
3736             BOOST_UBLAS_INLINE
operator -(const const_iterator2 & it) const3737             difference_type operator - (const const_iterator2 &it) const {
3738                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3739                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3740                 return it2_ - it.it2_;
3741             }
3742 
3743             // Dereference
3744             BOOST_UBLAS_INLINE
operator *() const3745             const_reference operator * () const {
3746                 // FIXME replace find with at_element
3747                 return (*this) ().data_ (*it1_, *it2_);
3748             }
3749 
3750 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
3751             BOOST_UBLAS_INLINE
3752 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3753             typename self_type::
3754 #endif
begin() const3755             const_iterator1 begin () const {
3756                 return const_iterator1 ((*this) (), it1_ ().begin (), it2_);
3757             }
3758             BOOST_UBLAS_INLINE
3759 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3760             typename self_type::
3761 #endif
end() const3762             const_iterator1 end () const {
3763                 return const_iterator1 ((*this) (), it1_ ().end (), it2_);
3764             }
3765             BOOST_UBLAS_INLINE
3766 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3767             typename self_type::
3768 #endif
rbegin() const3769             const_reverse_iterator1 rbegin () const {
3770                 return const_reverse_iterator1 (end ());
3771             }
3772             BOOST_UBLAS_INLINE
3773 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3774             typename self_type::
3775 #endif
rend() const3776             const_reverse_iterator1 rend () const {
3777                 return const_reverse_iterator1 (begin ());
3778             }
3779 #endif
3780 
3781             // Indices
3782             BOOST_UBLAS_INLINE
index1() const3783             size_type index1 () const {
3784                 return it1_.index ();
3785             }
3786             BOOST_UBLAS_INLINE
index2() const3787             size_type index2 () const {
3788                 return it2_.index ();
3789             }
3790 
3791             // Assignment
3792             BOOST_UBLAS_INLINE
operator =(const const_iterator2 & it)3793             const_iterator2 &operator = (const const_iterator2 &it) {
3794                 container_const_reference<self_type>::assign (&it ());
3795                 it1_ = it.it1_;
3796                 it2_ = it.it2_;
3797                 return *this;
3798             }
3799 
3800             // Comparison
3801             BOOST_UBLAS_INLINE
operator ==(const const_iterator2 & it) const3802             bool operator == (const const_iterator2 &it) const {
3803                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3804                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3805                 return it2_ == it.it2_;
3806             }
3807             BOOST_UBLAS_INLINE
operator <(const const_iterator2 & it) const3808             bool operator < (const const_iterator2 &it) const {
3809                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3810                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3811                 return it2_ < it.it2_;
3812             }
3813 
3814         private:
3815             const_subiterator1_type it1_;
3816             const_subiterator2_type it2_;
3817         };
3818 #endif
3819 
3820         BOOST_UBLAS_INLINE
begin2() const3821         const_iterator2 begin2 () const {
3822             return find2 (0, 0, 0);
3823         }
3824         BOOST_UBLAS_INLINE
end2() const3825         const_iterator2 end2 () const {
3826             return find2 (0, 0, size2 ());
3827         }
3828 
3829 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
3830         class iterator2:
3831             public container_reference<matrix_slice>,
3832             public iterator_base_traits<typename M::iterator2::iterator_category>::template
3833                         iterator_base<iterator2, value_type>::type {
3834         public:
3835             typedef typename M::iterator2::value_type value_type;
3836             typedef typename M::iterator2::difference_type difference_type;
3837             typedef typename M::reference reference;    //FIXME due to indexing access
3838             typedef typename M::iterator2::pointer pointer;
3839             typedef iterator1 dual_iterator_type;
3840             typedef reverse_iterator1 dual_reverse_iterator_type;
3841 
3842             // Construction and destruction
3843             BOOST_UBLAS_INLINE
iterator2()3844             iterator2 ():
3845                 container_reference<self_type> (), it1_ (), it2_ () {}
3846             BOOST_UBLAS_INLINE
iterator2(self_type & ms,const subiterator1_type & it1,const subiterator2_type & it2)3847             iterator2 (self_type &ms, const subiterator1_type &it1, const subiterator2_type &it2):
3848                 container_reference<self_type> (ms), it1_ (it1), it2_ (it2) {}
3849 
3850             // Arithmetic
3851             BOOST_UBLAS_INLINE
operator ++()3852             iterator2 &operator ++ () {
3853                 ++ it2_;
3854                 return *this;
3855             }
3856             BOOST_UBLAS_INLINE
operator --()3857             iterator2 &operator -- () {
3858                 -- it2_;
3859                 return *this;
3860             }
3861             BOOST_UBLAS_INLINE
operator +=(difference_type n)3862             iterator2 &operator += (difference_type n) {
3863                 it2_ += n;
3864                 return *this;
3865             }
3866             BOOST_UBLAS_INLINE
operator -=(difference_type n)3867             iterator2 &operator -= (difference_type n) {
3868                 it2_ -= n;
3869                 return *this;
3870             }
3871             BOOST_UBLAS_INLINE
operator -(const iterator2 & it) const3872             difference_type operator - (const iterator2 &it) const {
3873                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3874                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3875                 return it2_ - it.it2_;
3876             }
3877 
3878             // Dereference
3879             BOOST_UBLAS_INLINE
operator *() const3880             reference operator * () const {
3881                 // FIXME replace find with at_element
3882                 return (*this) ().data_ (*it1_, *it2_);
3883             }
3884 
3885 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
3886             BOOST_UBLAS_INLINE
3887 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3888             typename self_type::
3889 #endif
begin() const3890             iterator1 begin () const {
3891                 return iterator1 ((*this) (), it1_ ().begin (), it2_);
3892             }
3893             BOOST_UBLAS_INLINE
3894 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3895             typename self_type::
3896 #endif
end() const3897             iterator1 end () const {
3898                 return iterator1 ((*this) (), it1_ ().end (), it2_);
3899             }
3900             BOOST_UBLAS_INLINE
3901 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3902             typename self_type::
3903 #endif
rbegin() const3904             reverse_iterator1 rbegin () const {
3905                 return reverse_iterator1 (end ());
3906             }
3907             BOOST_UBLAS_INLINE
3908 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
3909             typename self_type::
3910 #endif
rend() const3911             reverse_iterator1 rend () const {
3912                 return reverse_iterator1 (begin ());
3913             }
3914 #endif
3915 
3916             // Indices
3917             BOOST_UBLAS_INLINE
index1() const3918             size_type index1 () const {
3919                 return it1_.index ();
3920             }
3921             BOOST_UBLAS_INLINE
index2() const3922             size_type index2 () const {
3923                 return it2_.index ();
3924             }
3925 
3926             // Assignment
3927             BOOST_UBLAS_INLINE
operator =(const iterator2 & it)3928             iterator2 &operator = (const iterator2 &it) {
3929                 container_reference<self_type>::assign (&it ());
3930                 it1_ = it.it1_;
3931                 it2_ = it.it2_;
3932                 return *this;
3933             }
3934 
3935             // Comparison
3936             BOOST_UBLAS_INLINE
operator ==(const iterator2 & it) const3937             bool operator == (const iterator2 &it) const {
3938                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3939                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3940                 return it2_ == it.it2_;
3941             }
3942             BOOST_UBLAS_INLINE
operator <(const iterator2 & it) const3943             bool operator < (const iterator2 &it) const {
3944                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
3945                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
3946                 return it2_ < it.it2_;
3947             }
3948 
3949         private:
3950             subiterator1_type it1_;
3951             subiterator2_type it2_;
3952 
3953             friend class const_iterator2;
3954         };
3955 #endif
3956 
3957         BOOST_UBLAS_INLINE
begin2()3958         iterator2 begin2 () {
3959             return find2 (0, 0, 0);
3960         }
3961         BOOST_UBLAS_INLINE
end2()3962         iterator2 end2 () {
3963             return find2 (0, 0, size2 ());
3964         }
3965 
3966         // Reverse iterators
3967 
3968         BOOST_UBLAS_INLINE
rbegin1() const3969         const_reverse_iterator1 rbegin1 () const {
3970             return const_reverse_iterator1 (end1 ());
3971         }
3972         BOOST_UBLAS_INLINE
rend1() const3973         const_reverse_iterator1 rend1 () const {
3974             return const_reverse_iterator1 (begin1 ());
3975         }
3976 
3977         BOOST_UBLAS_INLINE
rbegin1()3978         reverse_iterator1 rbegin1 () {
3979             return reverse_iterator1 (end1 ());
3980         }
3981         BOOST_UBLAS_INLINE
rend1()3982         reverse_iterator1 rend1 () {
3983             return reverse_iterator1 (begin1 ());
3984         }
3985 
3986         BOOST_UBLAS_INLINE
rbegin2() const3987         const_reverse_iterator2 rbegin2 () const {
3988             return const_reverse_iterator2 (end2 ());
3989         }
3990         BOOST_UBLAS_INLINE
rend2() const3991         const_reverse_iterator2 rend2 () const {
3992             return const_reverse_iterator2 (begin2 ());
3993         }
3994 
3995         BOOST_UBLAS_INLINE
rbegin2()3996         reverse_iterator2 rbegin2 () {
3997             return reverse_iterator2 (end2 ());
3998         }
3999         BOOST_UBLAS_INLINE
rend2()4000         reverse_iterator2 rend2 () {
4001             return reverse_iterator2 (begin2 ());
4002         }
4003 
4004     private:
4005         matrix_closure_type data_;
4006         slice_type s1_;
4007         slice_type s2_;
4008     };
4009 
4010     // Simple Projections
4011     template<class M>
4012     BOOST_UBLAS_INLINE
subslice(M & data,typename M::size_type start1,typename M::difference_type stride1,typename M::size_type size1,typename M::size_type start2,typename M::difference_type stride2,typename M::size_type size2)4013     matrix_slice<M> subslice (M &data, typename M::size_type start1, typename M::difference_type stride1, typename M::size_type size1, typename M::size_type start2, typename M::difference_type stride2, typename M::size_type size2) {
4014         typedef basic_slice<typename M::size_type, typename M::difference_type> slice_type;
4015         return matrix_slice<M> (data, slice_type (start1, stride1, size1), slice_type (start2, stride2, size2));
4016     }
4017     template<class M>
4018     BOOST_UBLAS_INLINE
subslice(const M & data,typename M::size_type start1,typename M::difference_type stride1,typename M::size_type size1,typename M::size_type start2,typename M::difference_type stride2,typename M::size_type size2)4019     matrix_slice<const M> subslice (const M &data, typename M::size_type start1, typename M::difference_type stride1, typename M::size_type size1, typename M::size_type start2, typename M::difference_type stride2, typename M::size_type size2) {
4020         typedef basic_slice<typename M::size_type, typename M::difference_type> slice_type;
4021         return matrix_slice<const M> (data (), slice_type (start1, stride1, size1), slice_type (start2, stride2, size2));
4022     }
4023 
4024     // Generic Projections
4025     template<class M>
4026     BOOST_UBLAS_INLINE
project(M & data,const typename matrix_slice<M>::slice_type & s1,const typename matrix_slice<M>::slice_type & s2)4027     matrix_slice<M> project (M &data, const typename matrix_slice<M>::slice_type &s1, const typename matrix_slice<M>::slice_type &s2) {
4028         return matrix_slice<M> (data, s1, s2);
4029     }
4030     template<class M>
4031     BOOST_UBLAS_INLINE
project(const M & data,const typename matrix_slice<M>::slice_type & s1,const typename matrix_slice<M>::slice_type & s2)4032     const matrix_slice<const M> project (const M &data, const typename matrix_slice<M>::slice_type &s1, const typename matrix_slice<M>::slice_type &s2) {
4033         // ISSUE was: return matrix_slice<M> (const_cast<M &> (data), s1, s2);
4034         return matrix_slice<const M> (data, s1, s2);
4035     }
4036     // ISSUE in the following two functions it would be logical to use matrix_slice<V>::range_type but this confuses VC7.1 and 8.0
4037     template<class M>
4038     BOOST_UBLAS_INLINE
project(matrix_slice<M> & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)4039     matrix_slice<M> project (matrix_slice<M> &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
4040         return data.project (r1, r2);
4041     }
4042     template<class M>
4043     BOOST_UBLAS_INLINE
project(const matrix_slice<M> & data,const typename matrix_range<M>::range_type & r1,const typename matrix_range<M>::range_type & r2)4044     const matrix_slice<M> project (const matrix_slice<M> &data, const typename matrix_range<M>::range_type &r1, const typename matrix_range<M>::range_type &r2) {
4045         return data.project (r1, r2);
4046     }
4047     template<class M>
4048     BOOST_UBLAS_INLINE
project(matrix_slice<M> & data,const typename matrix_slice<M>::slice_type & s1,const typename matrix_slice<M>::slice_type & s2)4049     matrix_slice<M> project (matrix_slice<M> &data, const typename matrix_slice<M>::slice_type &s1, const typename matrix_slice<M>::slice_type &s2) {
4050         return data.project (s1, s2);
4051     }
4052     template<class M>
4053     BOOST_UBLAS_INLINE
project(const matrix_slice<M> & data,const typename matrix_slice<M>::slice_type & s1,const typename matrix_slice<M>::slice_type & s2)4054     const matrix_slice<M> project (const matrix_slice<M> &data, const typename matrix_slice<M>::slice_type &s1, const typename matrix_slice<M>::slice_type &s2) {
4055         return data.project (s1, s2);
4056     }
4057 
4058     // Specialization of temporary_traits
4059     template <class M>
4060     struct matrix_temporary_traits< matrix_slice<M> >
4061     : matrix_temporary_traits< M > {};
4062 
4063     template <class M>
4064     struct vector_temporary_traits< matrix_slice<M> >
4065     : vector_temporary_traits< M > {};
4066 
4067     // Matrix based indirection class
4068     // Contributed by Toon Knapen.
4069     // Extended and optimized by Kresimir Fresl.
4070     template<class M, class IA>
4071     class matrix_indirect:
4072         public matrix_expression<matrix_indirect<M, IA> > {
4073 
4074         typedef matrix_indirect<M, IA> self_type;
4075     public:
4076 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
4077         using matrix_expression<self_type>::operator ();
4078 #endif
4079         typedef M matrix_type;
4080         typedef IA indirect_array_type;
4081         typedef typename M::size_type size_type;
4082         typedef typename M::difference_type difference_type;
4083         typedef typename M::value_type value_type;
4084         typedef typename M::const_reference const_reference;
4085         typedef typename boost::mpl::if_<boost::is_const<M>,
4086                                           typename M::const_reference,
4087                                           typename M::reference>::type reference;
4088         typedef typename boost::mpl::if_<boost::is_const<M>,
4089                                           typename M::const_closure_type,
4090                                           typename M::closure_type>::type matrix_closure_type;
4091         typedef basic_range<size_type, difference_type> range_type;
4092         typedef basic_slice<size_type, difference_type> slice_type;
4093         typedef const self_type const_closure_type;
4094         typedef self_type closure_type;
4095         typedef typename storage_restrict_traits<typename M::storage_category,
4096                                                  dense_proxy_tag>::storage_category storage_category;
4097         typedef typename M::orientation_category orientation_category;
4098 
4099         // Construction and destruction
4100         BOOST_UBLAS_INLINE
matrix_indirect(matrix_type & data,size_type size1,size_type size2)4101         matrix_indirect (matrix_type &data, size_type size1, size_type size2):
4102             data_ (data), ia1_ (size1), ia2_ (size2) {}
4103         BOOST_UBLAS_INLINE
matrix_indirect(matrix_type & data,const indirect_array_type & ia1,const indirect_array_type & ia2)4104         matrix_indirect (matrix_type &data, const indirect_array_type &ia1, const indirect_array_type &ia2):
4105             data_ (data), ia1_ (ia1.preprocess (data.size1 ())), ia2_ (ia2.preprocess (data.size2 ())) {}
4106         BOOST_UBLAS_INLINE
matrix_indirect(const matrix_closure_type & data,const indirect_array_type & ia1,const indirect_array_type & ia2,int)4107         matrix_indirect (const matrix_closure_type &data, const indirect_array_type &ia1, const indirect_array_type &ia2, int):
4108             data_ (data), ia1_ (ia1.preprocess (data.size1 ())), ia2_ (ia2.preprocess (data.size2 ())) {}
4109 
4110         // Accessors
4111         BOOST_UBLAS_INLINE
size1() const4112         size_type size1 () const {
4113             return ia1_.size ();
4114         }
4115         BOOST_UBLAS_INLINE
size2() const4116         size_type size2 () const {
4117             return ia2_.size ();
4118         }
4119         BOOST_UBLAS_INLINE
indirect1() const4120         const indirect_array_type &indirect1 () const {
4121             return ia1_;
4122         }
4123         BOOST_UBLAS_INLINE
indirect1()4124         indirect_array_type &indirect1 () {
4125             return ia1_;
4126         }
4127         BOOST_UBLAS_INLINE
indirect2() const4128         const indirect_array_type &indirect2 () const {
4129             return ia2_;
4130         }
4131         BOOST_UBLAS_INLINE
indirect2()4132         indirect_array_type &indirect2 () {
4133             return ia2_;
4134         }
4135 
4136         // Storage accessors
4137         BOOST_UBLAS_INLINE
data() const4138         const matrix_closure_type &data () const {
4139             return data_;
4140         }
4141         BOOST_UBLAS_INLINE
data()4142         matrix_closure_type &data () {
4143             return data_;
4144         }
4145 
4146         // Element access
4147 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
4148         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const4149         const_reference operator () (size_type i, size_type j) const {
4150             return data_ (ia1_ (i), ia2_ (j));
4151         }
4152         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j)4153         reference operator () (size_type i, size_type j) {
4154             return data_ (ia1_ (i), ia2_ (j));
4155         }
4156 #else
4157         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const4158         reference operator () (size_type i, size_type j) const {
4159             return data_ (ia1_ (i), ia2_ (j));
4160         }
4161 #endif
4162 
4163         // ISSUE can this be done in free project function?
4164         // Although a const function can create a non-const proxy to a non-const object
4165         // Critical is that matrix_type and data_ (vector_closure_type) are const correct
4166         BOOST_UBLAS_INLINE
project(const range_type & r1,const range_type & r2) const4167         matrix_indirect<matrix_type, indirect_array_type> project (const range_type &r1, const range_type &r2) const {
4168             return matrix_indirect<matrix_type, indirect_array_type> (data_, ia1_.compose (r1.preprocess (data_.size1 ())), ia2_.compose (r2.preprocess (data_.size2 ())), 0);
4169         }
4170         BOOST_UBLAS_INLINE
project(const slice_type & s1,const slice_type & s2) const4171         matrix_indirect<matrix_type, indirect_array_type> project (const slice_type &s1, const slice_type &s2) const {
4172             return matrix_indirect<matrix_type, indirect_array_type> (data_, ia1_.compose (s1.preprocess (data_.size1 ())), ia2_.compose (s2.preprocess (data_.size2 ())), 0);
4173         }
4174         BOOST_UBLAS_INLINE
project(const indirect_array_type & ia1,const indirect_array_type & ia2) const4175         matrix_indirect<matrix_type, indirect_array_type> project (const indirect_array_type &ia1, const indirect_array_type &ia2) const {
4176             return matrix_indirect<matrix_type, indirect_array_type> (data_, ia1_.compose (ia1.preprocess (data_.size1 ())), ia2_.compose (ia2.preprocess (data_.size2 ())), 0);
4177         }
4178 
4179         // Assignment
4180         BOOST_UBLAS_INLINE
operator =(const matrix_indirect & mi)4181         matrix_indirect &operator = (const matrix_indirect &mi) {
4182             matrix_assign<scalar_assign> (*this, mi);
4183             return *this;
4184         }
4185         BOOST_UBLAS_INLINE
assign_temporary(matrix_indirect & mi)4186         matrix_indirect &assign_temporary (matrix_indirect &mi) {
4187             return *this = mi;
4188         }
4189         template<class AE>
4190         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)4191         matrix_indirect &operator = (const matrix_expression<AE> &ae) {
4192             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (ae));
4193             return *this;
4194         }
4195         template<class AE>
4196         BOOST_UBLAS_INLINE
assign(const matrix_expression<AE> & ae)4197         matrix_indirect &assign (const matrix_expression<AE> &ae) {
4198             matrix_assign<scalar_assign> (*this, ae);
4199             return *this;
4200         }
4201         template<class AE>
4202         BOOST_UBLAS_INLINE
operator +=(const matrix_expression<AE> & ae)4203         matrix_indirect& operator += (const matrix_expression<AE> &ae) {
4204             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this + ae));
4205             return *this;
4206         }
4207         template<class AE>
4208         BOOST_UBLAS_INLINE
plus_assign(const matrix_expression<AE> & ae)4209         matrix_indirect &plus_assign (const matrix_expression<AE> &ae) {
4210             matrix_assign<scalar_plus_assign> (*this, ae);
4211             return *this;
4212         }
4213         template<class AE>
4214         BOOST_UBLAS_INLINE
operator -=(const matrix_expression<AE> & ae)4215         matrix_indirect& operator -= (const matrix_expression<AE> &ae) {
4216             matrix_assign<scalar_assign> (*this, typename matrix_temporary_traits<M>::type (*this - ae));
4217             return *this;
4218         }
4219         template<class AE>
4220         BOOST_UBLAS_INLINE
minus_assign(const matrix_expression<AE> & ae)4221         matrix_indirect &minus_assign (const matrix_expression<AE> &ae) {
4222             matrix_assign<scalar_minus_assign> (*this, ae);
4223             return *this;
4224         }
4225         template<class AT>
4226         BOOST_UBLAS_INLINE
operator *=(const AT & at)4227         matrix_indirect& operator *= (const AT &at) {
4228             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
4229             return *this;
4230         }
4231         template<class AT>
4232         BOOST_UBLAS_INLINE
operator /=(const AT & at)4233         matrix_indirect& operator /= (const AT &at) {
4234             matrix_assign_scalar<scalar_divides_assign> (*this, at);
4235             return *this;
4236         }
4237 
4238         // Closure comparison
4239         BOOST_UBLAS_INLINE
same_closure(const matrix_indirect & mi) const4240         bool same_closure (const matrix_indirect &mi) const {
4241             return (*this).data_.same_closure (mi.data_);
4242         }
4243 
4244         // Comparison
4245         BOOST_UBLAS_INLINE
operator ==(const matrix_indirect & mi) const4246         bool operator == (const matrix_indirect &mi) const {
4247             return (*this).data_ == mi.data_ && ia1_ == mi.ia1_ && ia2_ == mi.ia2_;
4248         }
4249 
4250         // Swapping
4251         BOOST_UBLAS_INLINE
swap(matrix_indirect mi)4252         void swap (matrix_indirect mi) {
4253             if (this != &mi) {
4254                 BOOST_UBLAS_CHECK (size1 () == mi.size1 (), bad_size ());
4255                 BOOST_UBLAS_CHECK (size2 () == mi.size2 (), bad_size ());
4256                 matrix_swap<scalar_swap> (*this, mi);
4257             }
4258         }
4259         BOOST_UBLAS_INLINE
swap(matrix_indirect mi1,matrix_indirect mi2)4260         friend void swap (matrix_indirect mi1, matrix_indirect mi2) {
4261             mi1.swap (mi2);
4262         }
4263 
4264         // Iterator types
4265     private:
4266         typedef typename IA::const_iterator const_subiterator1_type;
4267         typedef typename IA::const_iterator subiterator1_type;
4268         typedef typename IA::const_iterator const_subiterator2_type;
4269         typedef typename IA::const_iterator subiterator2_type;
4270 
4271     public:
4272 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
4273         typedef indexed_iterator1<matrix_indirect<matrix_type, indirect_array_type>,
4274                                   typename matrix_type::iterator1::iterator_category> iterator1;
4275         typedef indexed_iterator2<matrix_indirect<matrix_type, indirect_array_type>,
4276                                   typename matrix_type::iterator2::iterator_category> iterator2;
4277         typedef indexed_const_iterator1<matrix_indirect<matrix_type, indirect_array_type>,
4278                                         typename matrix_type::const_iterator1::iterator_category> const_iterator1;
4279         typedef indexed_const_iterator2<matrix_indirect<matrix_type, indirect_array_type>,
4280                                         typename matrix_type::const_iterator2::iterator_category> const_iterator2;
4281 #else
4282         class const_iterator1;
4283         class iterator1;
4284         class const_iterator2;
4285         class iterator2;
4286 #endif
4287         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
4288         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
4289         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
4290         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
4291 
4292         // Element lookup
4293         BOOST_UBLAS_INLINE
find1(int,size_type i,size_type j) const4294         const_iterator1 find1 (int /* rank */, size_type i, size_type j) const {
4295 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
4296             return const_iterator1 (*this, i, j);
4297 #else
4298             return const_iterator1 (*this, ia1_.begin () + i, ia2_.begin () + j);
4299 #endif
4300         }
4301         BOOST_UBLAS_INLINE
find1(int,size_type i,size_type j)4302         iterator1 find1 (int /* rank */, size_type i, size_type j) {
4303 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
4304             return iterator1 (*this, i, j);
4305 #else
4306             return iterator1 (*this, ia1_.begin () + i, ia2_.begin () + j);
4307 #endif
4308         }
4309         BOOST_UBLAS_INLINE
find2(int,size_type i,size_type j) const4310         const_iterator2 find2 (int /* rank */, size_type i, size_type j) const {
4311 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
4312             return const_iterator2 (*this, i, j);
4313 #else
4314             return const_iterator2 (*this, ia1_.begin () + i, ia2_.begin () + j);
4315 #endif
4316         }
4317         BOOST_UBLAS_INLINE
find2(int,size_type i,size_type j)4318         iterator2 find2 (int /* rank */, size_type i, size_type j) {
4319 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
4320             return iterator2 (*this, i, j);
4321 #else
4322             return iterator2 (*this, ia1_.begin () + i, ia2_.begin () + j);
4323 #endif
4324         }
4325 
4326         // Iterators simply are indices.
4327 
4328 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
4329         class const_iterator1:
4330             public container_const_reference<matrix_indirect>,
4331             public iterator_base_traits<typename M::const_iterator1::iterator_category>::template
4332                         iterator_base<const_iterator1, value_type>::type {
4333         public:
4334             typedef typename M::const_iterator1::value_type value_type;
4335             typedef typename M::const_iterator1::difference_type difference_type;
4336             typedef typename M::const_reference reference;    //FIXME due to indexing access
4337             typedef typename M::const_iterator1::pointer pointer;
4338             typedef const_iterator2 dual_iterator_type;
4339             typedef const_reverse_iterator2 dual_reverse_iterator_type;
4340 
4341             // Construction and destruction
4342             BOOST_UBLAS_INLINE
const_iterator1()4343             const_iterator1 ():
4344                 container_const_reference<self_type> (), it1_ (), it2_ () {}
4345             BOOST_UBLAS_INLINE
const_iterator1(const self_type & mi,const const_subiterator1_type & it1,const const_subiterator2_type & it2)4346             const_iterator1 (const self_type &mi, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
4347                 container_const_reference<self_type> (mi), it1_ (it1), it2_ (it2) {}
4348             BOOST_UBLAS_INLINE
const_iterator1(const iterator1 & it)4349             const_iterator1 (const iterator1 &it):
4350                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
4351 
4352             // Arithmetic
4353             BOOST_UBLAS_INLINE
operator ++()4354             const_iterator1 &operator ++ () {
4355                 ++ it1_;
4356                 return *this;
4357             }
4358             BOOST_UBLAS_INLINE
operator --()4359             const_iterator1 &operator -- () {
4360                 -- it1_;
4361                 return *this;
4362             }
4363             BOOST_UBLAS_INLINE
operator +=(difference_type n)4364             const_iterator1 &operator += (difference_type n) {
4365                 it1_ += n;
4366                 return *this;
4367             }
4368             BOOST_UBLAS_INLINE
operator -=(difference_type n)4369             const_iterator1 &operator -= (difference_type n) {
4370                 it1_ -= n;
4371                 return *this;
4372             }
4373             BOOST_UBLAS_INLINE
operator -(const const_iterator1 & it) const4374             difference_type operator - (const const_iterator1 &it) const {
4375                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4376                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4377                 return it1_ - it.it1_;
4378             }
4379 
4380             // Dereference
4381             BOOST_UBLAS_INLINE
operator *() const4382             const_reference operator * () const {
4383                 // FIXME replace find with at_element
4384                 return (*this) ().data_ (*it1_, *it2_);
4385             }
4386 
4387 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
4388             BOOST_UBLAS_INLINE
4389 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4390             typename self_type::
4391 #endif
begin() const4392             const_iterator2 begin () const {
4393                 return const_iterator2 ((*this) (), it1_, it2_ ().begin ());
4394             }
4395             BOOST_UBLAS_INLINE
4396 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4397             typename self_type::
4398 #endif
end() const4399             const_iterator2 end () const {
4400                 return const_iterator2 ((*this) (), it1_, it2_ ().end ());
4401             }
4402             BOOST_UBLAS_INLINE
4403 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4404             typename self_type::
4405 #endif
rbegin() const4406             const_reverse_iterator2 rbegin () const {
4407                 return const_reverse_iterator2 (end ());
4408             }
4409             BOOST_UBLAS_INLINE
4410 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4411             typename self_type::
4412 #endif
rend() const4413             const_reverse_iterator2 rend () const {
4414                 return const_reverse_iterator2 (begin ());
4415             }
4416 #endif
4417 
4418             // Indices
4419             BOOST_UBLAS_INLINE
index1() const4420             size_type index1 () const {
4421                 return it1_.index ();
4422             }
4423             BOOST_UBLAS_INLINE
index2() const4424             size_type index2 () const {
4425                 return it2_.index ();
4426             }
4427 
4428             // Assignment
4429             BOOST_UBLAS_INLINE
operator =(const const_iterator1 & it)4430             const_iterator1 &operator = (const const_iterator1 &it) {
4431                 container_const_reference<self_type>::assign (&it ());
4432                 it1_ = it.it1_;
4433                 it2_ = it.it2_;
4434                 return *this;
4435             }
4436 
4437             // Comparison
4438             BOOST_UBLAS_INLINE
operator ==(const const_iterator1 & it) const4439             bool operator == (const const_iterator1 &it) const {
4440                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4441                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4442                 return it1_ == it.it1_;
4443             }
4444             BOOST_UBLAS_INLINE
operator <(const const_iterator1 & it) const4445             bool operator < (const const_iterator1 &it) const {
4446                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4447                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4448                 return it1_ < it.it1_;
4449             }
4450 
4451         private:
4452             const_subiterator1_type it1_;
4453             const_subiterator2_type it2_;
4454         };
4455 #endif
4456 
4457         BOOST_UBLAS_INLINE
begin1() const4458         const_iterator1 begin1 () const {
4459             return find1 (0, 0, 0);
4460         }
4461         BOOST_UBLAS_INLINE
end1() const4462         const_iterator1 end1 () const {
4463             return find1 (0, size1 (), 0);
4464         }
4465 
4466 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
4467         class iterator1:
4468             public container_reference<matrix_indirect>,
4469             public iterator_base_traits<typename M::iterator1::iterator_category>::template
4470                         iterator_base<iterator1, value_type>::type {
4471         public:
4472             typedef typename M::iterator1::value_type value_type;
4473             typedef typename M::iterator1::difference_type difference_type;
4474             typedef typename M::reference reference;    //FIXME due to indexing access
4475             typedef typename M::iterator1::pointer pointer;
4476             typedef iterator2 dual_iterator_type;
4477             typedef reverse_iterator2 dual_reverse_iterator_type;
4478 
4479             // Construction and destruction
4480             BOOST_UBLAS_INLINE
iterator1()4481             iterator1 ():
4482                 container_reference<self_type> (), it1_ (), it2_ () {}
4483             BOOST_UBLAS_INLINE
iterator1(self_type & mi,const subiterator1_type & it1,const subiterator2_type & it2)4484             iterator1 (self_type &mi, const subiterator1_type &it1, const subiterator2_type &it2):
4485                 container_reference<self_type> (mi), it1_ (it1), it2_ (it2) {}
4486 
4487             // Arithmetic
4488             BOOST_UBLAS_INLINE
operator ++()4489             iterator1 &operator ++ () {
4490                 ++ it1_;
4491                 return *this;
4492             }
4493             BOOST_UBLAS_INLINE
operator --()4494             iterator1 &operator -- () {
4495                 -- it1_;
4496                 return *this;
4497             }
4498             BOOST_UBLAS_INLINE
operator +=(difference_type n)4499             iterator1 &operator += (difference_type n) {
4500                 it1_ += n;
4501                 return *this;
4502             }
4503             BOOST_UBLAS_INLINE
operator -=(difference_type n)4504             iterator1 &operator -= (difference_type n) {
4505                 it1_ -= n;
4506                 return *this;
4507             }
4508             BOOST_UBLAS_INLINE
operator -(const iterator1 & it) const4509             difference_type operator - (const iterator1 &it) const {
4510                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4511                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4512                 return it1_ - it.it1_;
4513             }
4514 
4515             // Dereference
4516             BOOST_UBLAS_INLINE
operator *() const4517             reference operator * () const {
4518                 // FIXME replace find with at_element
4519                 return (*this) ().data_ (*it1_, *it2_);
4520             }
4521 
4522 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
4523             BOOST_UBLAS_INLINE
4524 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4525             typename self_type::
4526 #endif
begin() const4527             iterator2 begin () const {
4528                 return iterator2 ((*this) (), it1_, it2_ ().begin ());
4529             }
4530             BOOST_UBLAS_INLINE
4531 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4532             typename self_type::
4533 #endif
end() const4534             iterator2 end () const {
4535                 return iterator2 ((*this) (), it1_, it2_ ().end ());
4536             }
4537             BOOST_UBLAS_INLINE
4538 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4539             typename self_type::
4540 #endif
rbegin() const4541             reverse_iterator2 rbegin () const {
4542                 return reverse_iterator2 (end ());
4543             }
4544             BOOST_UBLAS_INLINE
4545 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4546             typename self_type::
4547 #endif
rend() const4548             reverse_iterator2 rend () const {
4549                 return reverse_iterator2 (begin ());
4550             }
4551 #endif
4552 
4553             // Indices
4554             BOOST_UBLAS_INLINE
index1() const4555             size_type index1 () const {
4556                 return it1_.index ();
4557             }
4558             BOOST_UBLAS_INLINE
index2() const4559             size_type index2 () const {
4560                 return it2_.index ();
4561             }
4562 
4563             // Assignment
4564             BOOST_UBLAS_INLINE
operator =(const iterator1 & it)4565             iterator1 &operator = (const iterator1 &it) {
4566                 container_reference<self_type>::assign (&it ());
4567                 it1_ = it.it1_;
4568                 it2_ = it.it2_;
4569                 return *this;
4570             }
4571 
4572             // Comparison
4573             BOOST_UBLAS_INLINE
operator ==(const iterator1 & it) const4574             bool operator == (const iterator1 &it) const {
4575                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4576                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4577                 return it1_ == it.it1_;
4578             }
4579             BOOST_UBLAS_INLINE
operator <(const iterator1 & it) const4580             bool operator < (const iterator1 &it) const {
4581                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4582                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
4583                 return it1_ < it.it1_;
4584             }
4585 
4586         private:
4587             subiterator1_type it1_;
4588             subiterator2_type it2_;
4589 
4590             friend class const_iterator1;
4591         };
4592 #endif
4593 
4594         BOOST_UBLAS_INLINE
begin1()4595         iterator1 begin1 () {
4596             return find1 (0, 0, 0);
4597         }
4598         BOOST_UBLAS_INLINE
end1()4599         iterator1 end1 () {
4600             return find1 (0, size1 (), 0);
4601         }
4602 
4603 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
4604         class const_iterator2:
4605             public container_const_reference<matrix_indirect>,
4606             public iterator_base_traits<typename M::const_iterator2::iterator_category>::template
4607                         iterator_base<const_iterator2, value_type>::type {
4608         public:
4609             typedef typename M::const_iterator2::value_type value_type;
4610             typedef typename M::const_iterator2::difference_type difference_type;
4611             typedef typename M::const_reference reference;    //FIXME due to indexing access
4612             typedef typename M::const_iterator2::pointer pointer;
4613             typedef const_iterator1 dual_iterator_type;
4614             typedef const_reverse_iterator1 dual_reverse_iterator_type;
4615 
4616             // Construction and destruction
4617             BOOST_UBLAS_INLINE
const_iterator2()4618             const_iterator2 ():
4619                 container_const_reference<self_type> (), it1_ (), it2_ () {}
4620             BOOST_UBLAS_INLINE
const_iterator2(const self_type & mi,const const_subiterator1_type & it1,const const_subiterator2_type & it2)4621             const_iterator2 (const self_type &mi, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
4622                 container_const_reference<self_type> (mi), it1_ (it1), it2_ (it2) {}
4623             BOOST_UBLAS_INLINE
const_iterator2(const iterator2 & it)4624             const_iterator2 (const iterator2 &it):
4625                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
4626 
4627             // Arithmetic
4628             BOOST_UBLAS_INLINE
operator ++()4629             const_iterator2 &operator ++ () {
4630                 ++ it2_;
4631                 return *this;
4632             }
4633             BOOST_UBLAS_INLINE
operator --()4634             const_iterator2 &operator -- () {
4635                 -- it2_;
4636                 return *this;
4637             }
4638             BOOST_UBLAS_INLINE
operator +=(difference_type n)4639             const_iterator2 &operator += (difference_type n) {
4640                 it2_ += n;
4641                 return *this;
4642             }
4643             BOOST_UBLAS_INLINE
operator -=(difference_type n)4644             const_iterator2 &operator -= (difference_type n) {
4645                 it2_ -= n;
4646                 return *this;
4647             }
4648             BOOST_UBLAS_INLINE
operator -(const const_iterator2 & it) const4649             difference_type operator - (const const_iterator2 &it) const {
4650                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4651                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4652                 return it2_ - it.it2_;
4653             }
4654 
4655             // Dereference
4656             BOOST_UBLAS_INLINE
operator *() const4657             const_reference operator * () const {
4658                 // FIXME replace find with at_element
4659                 return (*this) ().data_ (*it1_, *it2_);
4660             }
4661 
4662 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
4663             BOOST_UBLAS_INLINE
4664 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4665             typename self_type::
4666 #endif
begin() const4667             const_iterator1 begin () const {
4668                 return const_iterator1 ((*this) (), it1_ ().begin (), it2_);
4669             }
4670             BOOST_UBLAS_INLINE
4671 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4672             typename self_type::
4673 #endif
end() const4674             const_iterator1 end () const {
4675                 return const_iterator1 ((*this) (), it1_ ().end (), it2_);
4676             }
4677             BOOST_UBLAS_INLINE
4678 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4679             typename self_type::
4680 #endif
rbegin() const4681             const_reverse_iterator1 rbegin () const {
4682                 return const_reverse_iterator1 (end ());
4683             }
4684             BOOST_UBLAS_INLINE
4685 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4686             typename self_type::
4687 #endif
rend() const4688             const_reverse_iterator1 rend () const {
4689                 return const_reverse_iterator1 (begin ());
4690             }
4691 #endif
4692 
4693             // Indices
4694             BOOST_UBLAS_INLINE
index1() const4695             size_type index1 () const {
4696                 return it1_.index ();
4697             }
4698             BOOST_UBLAS_INLINE
index2() const4699             size_type index2 () const {
4700                 return it2_.index ();
4701             }
4702 
4703             // Assignment
4704             BOOST_UBLAS_INLINE
operator =(const const_iterator2 & it)4705             const_iterator2 &operator = (const const_iterator2 &it) {
4706                 container_const_reference<self_type>::assign (&it ());
4707                 it1_ = it.it1_;
4708                 it2_ = it.it2_;
4709                 return *this;
4710             }
4711 
4712             // Comparison
4713             BOOST_UBLAS_INLINE
operator ==(const const_iterator2 & it) const4714             bool operator == (const const_iterator2 &it) const {
4715                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4716                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4717                 return it2_ == it.it2_;
4718             }
4719             BOOST_UBLAS_INLINE
operator <(const const_iterator2 & it) const4720             bool operator < (const const_iterator2 &it) const {
4721                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4722                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4723                 return it2_ < it.it2_;
4724             }
4725 
4726         private:
4727             const_subiterator1_type it1_;
4728             const_subiterator2_type it2_;
4729         };
4730 #endif
4731 
4732         BOOST_UBLAS_INLINE
begin2() const4733         const_iterator2 begin2 () const {
4734             return find2 (0, 0, 0);
4735         }
4736         BOOST_UBLAS_INLINE
end2() const4737         const_iterator2 end2 () const {
4738             return find2 (0, 0, size2 ());
4739         }
4740 
4741 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
4742         class iterator2:
4743             public container_reference<matrix_indirect>,
4744             public iterator_base_traits<typename M::iterator2::iterator_category>::template
4745                         iterator_base<iterator2, value_type>::type {
4746         public:
4747             typedef typename M::iterator2::value_type value_type;
4748             typedef typename M::iterator2::difference_type difference_type;
4749             typedef typename M::reference reference;    //FIXME due to indexing access
4750             typedef typename M::iterator2::pointer pointer;
4751             typedef iterator1 dual_iterator_type;
4752             typedef reverse_iterator1 dual_reverse_iterator_type;
4753 
4754             // Construction and destruction
4755             BOOST_UBLAS_INLINE
iterator2()4756             iterator2 ():
4757                 container_reference<self_type> (), it1_ (), it2_ () {}
4758             BOOST_UBLAS_INLINE
iterator2(self_type & mi,const subiterator1_type & it1,const subiterator2_type & it2)4759             iterator2 (self_type &mi, const subiterator1_type &it1, const subiterator2_type &it2):
4760                 container_reference<self_type> (mi), it1_ (it1), it2_ (it2) {}
4761 
4762             // Arithmetic
4763             BOOST_UBLAS_INLINE
operator ++()4764             iterator2 &operator ++ () {
4765                 ++ it2_;
4766                 return *this;
4767             }
4768             BOOST_UBLAS_INLINE
operator --()4769             iterator2 &operator -- () {
4770                 -- it2_;
4771                 return *this;
4772             }
4773             BOOST_UBLAS_INLINE
operator +=(difference_type n)4774             iterator2 &operator += (difference_type n) {
4775                 it2_ += n;
4776                 return *this;
4777             }
4778             BOOST_UBLAS_INLINE
operator -=(difference_type n)4779             iterator2 &operator -= (difference_type n) {
4780                 it2_ -= n;
4781                 return *this;
4782             }
4783             BOOST_UBLAS_INLINE
operator -(const iterator2 & it) const4784             difference_type operator - (const iterator2 &it) const {
4785                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4786                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4787                 return it2_ - it.it2_;
4788             }
4789 
4790             // Dereference
4791             BOOST_UBLAS_INLINE
operator *() const4792             reference operator * () const {
4793                 // FIXME replace find with at_element
4794                 return (*this) ().data_ (*it1_, *it2_);
4795             }
4796 
4797 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
4798             BOOST_UBLAS_INLINE
4799 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4800             typename self_type::
4801 #endif
begin() const4802             iterator1 begin () const {
4803                 return iterator1 ((*this) (), it1_ ().begin (), it2_);
4804             }
4805             BOOST_UBLAS_INLINE
4806 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4807             typename self_type::
4808 #endif
end() const4809             iterator1 end () const {
4810                 return iterator1 ((*this) (), it1_ ().end (), it2_);
4811             }
4812             BOOST_UBLAS_INLINE
4813 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4814             typename self_type::
4815 #endif
rbegin() const4816             reverse_iterator1 rbegin () const {
4817                 return reverse_iterator1 (end ());
4818             }
4819             BOOST_UBLAS_INLINE
4820 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
4821             typename self_type::
4822 #endif
rend() const4823             reverse_iterator1 rend () const {
4824                 return reverse_iterator1 (begin ());
4825             }
4826 #endif
4827 
4828             // Indices
4829             BOOST_UBLAS_INLINE
index1() const4830             size_type index1 () const {
4831                 return it1_.index ();
4832             }
4833             BOOST_UBLAS_INLINE
index2() const4834             size_type index2 () const {
4835                 return it2_.index ();
4836             }
4837 
4838             // Assignment
4839             BOOST_UBLAS_INLINE
operator =(const iterator2 & it)4840             iterator2 &operator = (const iterator2 &it) {
4841                 container_reference<self_type>::assign (&it ());
4842                 it1_ = it.it1_;
4843                 it2_ = it.it2_;
4844                 return *this;
4845             }
4846 
4847             // Comparison
4848             BOOST_UBLAS_INLINE
operator ==(const iterator2 & it) const4849             bool operator == (const iterator2 &it) const {
4850                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4851                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4852                 return it2_ == it.it2_;
4853             }
4854             BOOST_UBLAS_INLINE
operator <(const iterator2 & it) const4855             bool operator < (const iterator2 &it) const {
4856                 BOOST_UBLAS_CHECK ((*this) ().same_closure  (it ()), external_logic ());
4857                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
4858                 return it2_ < it.it2_;
4859             }
4860 
4861         private:
4862             subiterator1_type it1_;
4863             subiterator2_type it2_;
4864 
4865             friend class const_iterator2;
4866         };
4867 #endif
4868 
4869         BOOST_UBLAS_INLINE
begin2()4870         iterator2 begin2 () {
4871             return find2 (0, 0, 0);
4872         }
4873         BOOST_UBLAS_INLINE
end2()4874         iterator2 end2 () {
4875             return find2 (0, 0, size2 ());
4876         }
4877 
4878         // Reverse iterators
4879 
4880         BOOST_UBLAS_INLINE
rbegin1() const4881         const_reverse_iterator1 rbegin1 () const {
4882             return const_reverse_iterator1 (end1 ());
4883         }
4884         BOOST_UBLAS_INLINE
rend1() const4885         const_reverse_iterator1 rend1 () const {
4886             return const_reverse_iterator1 (begin1 ());
4887         }
4888 
4889         BOOST_UBLAS_INLINE
rbegin1()4890         reverse_iterator1 rbegin1 () {
4891             return reverse_iterator1 (end1 ());
4892         }
4893         BOOST_UBLAS_INLINE
rend1()4894         reverse_iterator1 rend1 () {
4895             return reverse_iterator1 (begin1 ());
4896         }
4897 
4898         BOOST_UBLAS_INLINE
rbegin2() const4899         const_reverse_iterator2 rbegin2 () const {
4900             return const_reverse_iterator2 (end2 ());
4901         }
4902         BOOST_UBLAS_INLINE
rend2() const4903         const_reverse_iterator2 rend2 () const {
4904             return const_reverse_iterator2 (begin2 ());
4905         }
4906 
4907         BOOST_UBLAS_INLINE
rbegin2()4908         reverse_iterator2 rbegin2 () {
4909             return reverse_iterator2 (end2 ());
4910         }
4911         BOOST_UBLAS_INLINE
rend2()4912         reverse_iterator2 rend2 () {
4913             return reverse_iterator2 (begin2 ());
4914         }
4915 
4916     private:
4917         matrix_closure_type data_;
4918         indirect_array_type ia1_;
4919         indirect_array_type ia2_;
4920     };
4921 
4922     // Projections
4923     template<class M, class A>
4924     BOOST_UBLAS_INLINE
project(M & data,const indirect_array<A> & ia1,const indirect_array<A> & ia2)4925     matrix_indirect<M, indirect_array<A> > project (M &data, const indirect_array<A> &ia1, const indirect_array<A> &ia2) {
4926         return matrix_indirect<M, indirect_array<A> > (data, ia1, ia2);
4927     }
4928     template<class M, class A>
4929     BOOST_UBLAS_INLINE
project(const M & data,const indirect_array<A> & ia1,const indirect_array<A> & ia2)4930     const matrix_indirect<const M, indirect_array<A> > project (const M &data, const indirect_array<A> &ia1, const indirect_array<A> &ia2) {
4931         // ISSUE was: return matrix_indirect<M, indirect_array<A> > (const_cast<M &> (data), ia1, ia2);
4932         return matrix_indirect<const M, indirect_array<A> > (data, ia1, ia2);
4933     }
4934     template<class M, class IA>
4935     BOOST_UBLAS_INLINE
project(matrix_indirect<M,IA> & data,const typename matrix_indirect<M,IA>::range_type & r1,const typename matrix_indirect<M,IA>::range_type & r2)4936     matrix_indirect<M, IA> project (matrix_indirect<M, IA> &data, const typename matrix_indirect<M, IA>::range_type &r1, const typename matrix_indirect<M, IA>::range_type &r2) {
4937         return data.project (r1, r2);
4938     }
4939     template<class M, class IA>
4940     BOOST_UBLAS_INLINE
project(const matrix_indirect<M,IA> & data,const typename matrix_indirect<M,IA>::range_type & r1,const typename matrix_indirect<M,IA>::range_type & r2)4941     const matrix_indirect<M, IA> project (const matrix_indirect<M, IA> &data, const typename matrix_indirect<M, IA>::range_type &r1, const typename matrix_indirect<M, IA>::range_type &r2) {
4942         return data.project (r1, r2);
4943     }
4944     template<class M, class IA>
4945     BOOST_UBLAS_INLINE
project(matrix_indirect<M,IA> & data,const typename matrix_indirect<M,IA>::slice_type & s1,const typename matrix_indirect<M,IA>::slice_type & s2)4946     matrix_indirect<M, IA> project (matrix_indirect<M, IA> &data, const typename matrix_indirect<M, IA>::slice_type &s1, const typename matrix_indirect<M, IA>::slice_type &s2) {
4947         return data.project (s1, s2);
4948     }
4949     template<class M, class IA>
4950     BOOST_UBLAS_INLINE
project(const matrix_indirect<M,IA> & data,const typename matrix_indirect<M,IA>::slice_type & s1,const typename matrix_indirect<M,IA>::slice_type & s2)4951     const matrix_indirect<M, IA> project (const matrix_indirect<M, IA> &data, const typename matrix_indirect<M, IA>::slice_type &s1, const typename matrix_indirect<M, IA>::slice_type &s2) {
4952         return data.project (s1, s2);
4953     }
4954     template<class M, class A>
4955     BOOST_UBLAS_INLINE
project(matrix_indirect<M,indirect_array<A>> & data,const indirect_array<A> & ia1,const indirect_array<A> & ia2)4956     matrix_indirect<M, indirect_array<A> > project (matrix_indirect<M, indirect_array<A> > &data, const indirect_array<A> &ia1, const indirect_array<A> &ia2) {
4957         return data.project (ia1, ia2);
4958     }
4959     template<class M, class A>
4960     BOOST_UBLAS_INLINE
project(const matrix_indirect<M,indirect_array<A>> & data,const indirect_array<A> & ia1,const indirect_array<A> & ia2)4961     const matrix_indirect<M, indirect_array<A> > project (const matrix_indirect<M, indirect_array<A> > &data, const indirect_array<A> &ia1, const indirect_array<A> &ia2) {
4962         return data.project (ia1, ia2);
4963     }
4964 
4965     /// Specialization of temporary_traits
4966     template <class M>
4967     struct matrix_temporary_traits< matrix_indirect<M> >
4968     : matrix_temporary_traits< M > {};
4969 
4970     template <class M>
4971     struct vector_temporary_traits< matrix_indirect<M> >
4972     : vector_temporary_traits< M > {};
4973 
4974 }}}
4975 
4976 #endif
4977