1 //
2 //  Copyright (c) 2000-2013
3 //  Joerg Walter, Mathias Koch, Athanasios Iliopoulos
4 //
5 //  Distributed under the Boost Software License, Version 1.0. (See
6 //  accompanying file LICENSE_1_0.txt or copy at
7 //  http://www.boost.org/LICENSE_1_0.txt)
8 //
9 //  The authors gratefully acknowledge the support of
10 //  GeNeSys mbH & Co. KG in producing this work.
11 //
12 
13 #ifndef _BOOST_UBLAS_BANDED_
14 #define _BOOST_UBLAS_BANDED_
15 
16 #include <boost/numeric/ublas/matrix.hpp>
17 #include <boost/numeric/ublas/detail/temporary.hpp>
18 
19 // Iterators based on ideas of Jeremy Siek
20 
21 namespace boost { namespace numeric { namespace ublas {
22 
23 
24 namespace hidden {
25 
26 
27 
28 /** \brief A helper for band_matrix indexing.
29  *
30  * The indexing happens as per the netlib description: http://www.netlib.org/lapack/lug/node124.html.
31  * In the case of a row_major matrix a different approach is followed;
32  */
33 template <class LayoutType>
34 class banded_indexing { };
35 
36 /** \brief A helper for indexing column major banded matrices.
37  *
38  */
39 template <>
40 class banded_indexing<column_major_tag> {
41 public:
42 
43     template <class T>
size(T,T size2)44     BOOST_UBLAS_INLINE static T size(T /*size1*/, T size2) {
45         return size2;
46     }
47 
48 //    template <class T>
49 //    BOOST_UBLAS_INLINE static bool valid_index(T size1, T /*size2*/, T lower, T upper, T i, T j) {
50 //        return (upper+i >= j) && i <= std::min(size1 - 1, j + lower); // upper + i is used by get_index. Maybe find a way to consolidate the operations to increase performance
51 //    }
52 
53     template <class T>
get_index(T,T size2,T lower,T upper,T i,T j)54     BOOST_UBLAS_INLINE static T get_index(T /*size1*/, T size2, T lower, T upper, T i, T j) {
55         return column_major::element (upper + i - j, lower + 1 + upper, j, size2);
56     }
57 };
58 
59 /** \brief A helper for indexing row major banded matrices.
60  *
61  */
62 template <>
63 class banded_indexing<row_major_tag> {
64 public:
65 
66     template <class T>
size(T size1,T)67     BOOST_UBLAS_INLINE static T size(T size1, T /*size2*/) {
68         return size1;
69     }
70 
71   //  template <class T>
72   //  BOOST_UBLAS_INLINE static bool valid_index(T /*size1*/, T  size2, T lower, T upper, T i, T j) {
73   //      return (lower+j >= i) && j <= std::min(size2 - 1, i + upper); // lower + j is used by get_index. Maybe find a way to consolidate the operations to increase performance
74   //  }
75 
76     template <class T>
get_index(T size1,T,T lower,T upper,T i,T j)77     BOOST_UBLAS_INLINE static T get_index(T size1, T /*size2*/, T lower, T upper, T i, T j) {
78         return row_major::element (i, size1, lower + j - i, lower + 1 + upper);
79     }
80 };
81 
82 }
83 
84     /** \brief A banded matrix of values of type \c T.
85      *
86      * For a \f$(mxn)\f$-dimensional banded matrix with \f$l\f$ lower and \f$u\f$ upper diagonals and
87      * \f$0 \leq i < m\f$ and \f$0 \leq j < n\f$, if \f$i>j+l\f$ or \f$i<j-u\f$ then \f$b_{i,j}=0\f$.
88      * The default storage for banded matrices is packed. Orientation and storage can also be specified.
89      * Default is \c row_major and and unbounded_array. It is \b not required by the storage to initialize
90      * elements of the matrix.
91      *
92      * \tparam T the type of object stored in the matrix (like double, float, complex, etc...)
93      * \tparam L the storage organization. It can be either \c row_major or \c column_major. Default is \c row_major
94      * \tparam A the type of Storage array. Default is \c unbounded_array
95      */
96     template<class T, class L, class A>
97     class banded_matrix:
98         public matrix_container<banded_matrix<T, L, A> > {
99 
100         typedef T *pointer;
101         typedef L layout_type;
102         typedef banded_matrix<T, L, A> self_type;
103 
104 
105 
106     public:
107 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
108         using matrix_container<self_type>::operator ();
109 #endif
110         typedef typename A::size_type size_type;
111         typedef typename A::difference_type difference_type;
112         typedef T value_type;
113         typedef const T &const_reference;
114         typedef T &reference;
115         typedef A array_type;
116         typedef const matrix_reference<const self_type> const_closure_type;
117         typedef matrix_reference<self_type> closure_type;
118         typedef vector<T, A> vector_temporary_type;
119         typedef matrix<T, L, A> matrix_temporary_type;  // general sub-matrix
120         typedef packed_tag storage_category;
121         typedef typename L::orientation_category orientation_category;
122 
123     private:
124     public:
125 
126         // Construction and destruction
127         BOOST_UBLAS_INLINE
banded_matrix()128         banded_matrix ():
129             matrix_container<self_type> (),
130             size1_ (0), size2_ (0),
131             lower_ (0), upper_ (0), data_ (0) {}
132         BOOST_UBLAS_INLINE
banded_matrix(size_type size1,size_type size2,size_type lower=0,size_type upper=0)133         banded_matrix (size_type size1, size_type size2, size_type lower = 0, size_type upper = 0):
134             matrix_container<self_type> (),
135             size1_ (size1), size2_ (size2),
136             lower_ (lower), upper_ (upper),
137 #if defined(BOOST_UBLAS_OWN_BANDED) || (BOOST_UBLAS_LEGACY_BANDED)
138             data_ ((std::max) (size1, size2) * (lower + 1 + upper))
139 #else
140             data_ ( hidden::banded_indexing<orientation_category>::size(size1, size2) * (lower + 1 + upper)) // This is the netlib layout as described here: http://www.netlib.org/lapack/lug/node124.html
141 #endif
142         {
143         }
144         BOOST_UBLAS_INLINE
banded_matrix(size_type size1,size_type size2,size_type lower,size_type upper,const array_type & data)145         banded_matrix (size_type size1, size_type size2, size_type lower, size_type upper, const array_type &data):
146             matrix_container<self_type> (),
147             size1_ (size1), size2_ (size2),
148             lower_ (lower), upper_ (upper), data_ (data) {}
149         BOOST_UBLAS_INLINE
banded_matrix(const banded_matrix & m)150         banded_matrix (const banded_matrix &m):
151             matrix_container<self_type> (),
152             size1_ (m.size1_), size2_ (m.size2_),
153             lower_ (m.lower_), upper_ (m.upper_), data_ (m.data_) {}
154         template<class AE>
155         BOOST_UBLAS_INLINE
banded_matrix(const matrix_expression<AE> & ae,size_type lower=0,size_type upper=0)156         banded_matrix (const matrix_expression<AE> &ae, size_type lower = 0, size_type upper = 0):
157             matrix_container<self_type> (),
158             size1_ (ae ().size1 ()), size2_ (ae ().size2 ()),
159             lower_ (lower), upper_ (upper),
160 #if defined(BOOST_UBLAS_OWN_BANDED) || (BOOST_UBLAS_LEGACY_BANDED)
161             data_ ((std::max) (size1_, size2_) * (lower_ + 1 + upper_))
162 #else
163             data_ ( hidden::banded_indexing<orientation_category>::size(size1_, size2_) * (lower_ + 1 + upper_)) // This is the netlib layout as described here: http://www.netlib.org/lapack/lug/node124.html
164 #endif
165         {
166             matrix_assign<scalar_assign> (*this, ae);
167         }
168 
169         // Accessors
170         BOOST_UBLAS_INLINE
size1() const171         size_type size1 () const {
172             return size1_;
173         }
174         BOOST_UBLAS_INLINE
size2() const175         size_type size2 () const {
176             return size2_;
177         }
178         BOOST_UBLAS_INLINE
lower() const179         size_type lower () const {
180             return lower_;
181         }
182         BOOST_UBLAS_INLINE
upper() const183         size_type upper () const {
184             return upper_;
185         }
186 
187         // Storage accessors
188         BOOST_UBLAS_INLINE
data() const189         const array_type &data () const {
190             return data_;
191         }
192         BOOST_UBLAS_INLINE
data()193         array_type &data () {
194             return data_;
195         }
196 
197 #if !defined (BOOST_UBLAS_OWN_BANDED)||(BOOST_UBLAS_LEGACY_BANDED)
198         BOOST_UBLAS_INLINE
is_element_in_band(size_type i,size_type j) const199         bool is_element_in_band(size_type i, size_type j) const{
200             //return (upper_+i >= j) && i <= std::min(size1() - 1, j + lower_); // We don't need to check if i is outside because it is checked anyway in the accessors.
201             return (upper_+i >= j) && i <= ( j + lower_); // Essentially this band has "infinite" positive dimensions
202         }
203 #endif
204         // Resizing
205         BOOST_UBLAS_INLINE
resize(size_type size1,size_type size2,size_type lower=0,size_type upper=0,bool preserve=true)206         void resize (size_type size1, size_type size2, size_type lower = 0, size_type upper = 0, bool preserve = true) {
207             if (preserve) {
208                 self_type temporary (size1, size2, lower, upper);
209                 detail::matrix_resize_preserve<layout_type> (*this, temporary);
210             }
211             else {
212                 data ().resize ((std::max) (size1, size2) * (lower + 1 + upper));
213                 size1_ = size1;
214                 size2_ = size2;
215                 lower_ = lower;
216                 upper_ = upper;
217             }
218         }
219 
220         BOOST_UBLAS_INLINE
resize_packed_preserve(size_type size1,size_type size2,size_type lower=0,size_type upper=0)221         void resize_packed_preserve (size_type size1, size_type size2, size_type lower = 0, size_type upper = 0) {
222             size1_ = size1;
223             size2_ = size2;
224             lower_ = lower;
225             upper_ = upper;
226             data ().resize ((std::max) (size1, size2) * (lower + 1 + upper), value_type ());
227         }
228 
229         // Element access
230         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const231         const_reference operator () (size_type i, size_type j) const {
232             BOOST_UBLAS_CHECK (i < size1_, bad_index ());
233             BOOST_UBLAS_CHECK (j < size2_, bad_index ());
234 #ifdef BOOST_UBLAS_OWN_BANDED
235             const size_type k = (std::max) (i, j);
236             const size_type l = lower_ + j - i;
237             if (k < (std::max) (size1_, size2_) && // TODO: probably use BOOST_UBLAS_CHECK here instead of if
238                 l < lower_ + 1 + upper_)
239                 return data () [layout_type::element (k, (std::max) (size1_, size2_),
240                                                        l, lower_ + 1 + upper_)];
241 #elif BOOST_UBLAS_LEGACY_BANDED // Prior to version: TODO: add version this is actually incorporated in
242             const size_type k = j;
243             const size_type l = upper_ + i - j;
244             if (k < size2_ &&
245                 l < lower_ + 1 + upper_)
246                 return data () [layout_type::element (k, size2_,
247                                                        l, lower_ + 1 + upper_)];
248 #else  // New default
249             // This is the netlib layout as described here: http://www.netlib.org/lapack/lug/node124.html
250             if ( is_element_in_band( i, j) ) {
251                 return data () [hidden::banded_indexing<orientation_category>::get_index(size1_, size2_, lower_, upper_, i, j)];
252             }
253 #endif
254             return zero_;
255         }
256 
257         BOOST_UBLAS_INLINE
at_element(size_type i,size_type j)258         reference at_element (size_type i, size_type j) {
259             BOOST_UBLAS_CHECK (i < size1_, bad_index ());
260             BOOST_UBLAS_CHECK (j < size2_, bad_index ());
261 #ifdef BOOST_UBLAS_OWN_BANDED
262             const size_type k = (std::max) (i, j);
263             const size_type l = lower_ + j - i; // TODO: Don't we need an if or BOOST_UBLAS_CHECK HERE?
264             return data () [layout_type::element (k, (std::max) (size1_, size2_),
265                                                    l, lower_ + 1 + upper_)];
266 #elif BOOST_UBLAS_LEGACY_BANDED // Prior to version: TODO: add version this is actually incorporated in
267             const size_type k = j;
268             const size_type l = upper_ + i - j;
269             if (! (k < size2_ &&
270                    l < lower_ + 1 + upper_) ) {
271                 bad_index ().raise ();
272                 // NEVER reached
273             }
274             return data () [layout_type::element (k, size2_,
275                                                        l, lower_ + 1 + upper_)];
276 #else
277             // This is the netlib layout as described here: http://www.netlib.org/lapack/lug/node124.html
278             BOOST_UBLAS_CHECK(is_element_in_band( i, j) , bad_index());
279             return data () [hidden::banded_indexing<orientation_category>::get_index(size1_, size2_, lower_, upper_, i, j)];
280 #endif
281         }
282         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j)283         reference operator () (size_type i, size_type j) {
284             BOOST_UBLAS_CHECK (i < size1_, bad_index ());
285             BOOST_UBLAS_CHECK (j < size2_, bad_index ());
286 #ifdef BOOST_UBLAS_OWN_BANDED
287             const size_type k = (std::max) (i, j);
288             const size_type l = lower_ + j - i;
289             if (! (k < (std::max) (size1_, size2_) && // TODO: probably use BOOST_UBLAS_CHECK here instead of if
290                   l < lower_ + 1 + upper_) ) {
291                 bad_index ().raise ();
292                 // NEVER reached
293             }
294             return data () [layout_type::element (k, (std::max) (size1_, size2_),
295                                                        l, lower_ + 1 + upper_)];
296 #elif BOOST_UBLAS_LEGACY_BANDED // Prior to version: TODO: add version this is actually incorporated in
297             const size_type k = j;
298             const size_type l = upper_ + i - j;
299             if (! (k < size2_ &&
300                    l < lower_ + 1 + upper_) ) {
301                 bad_index ().raise ();
302                 // NEVER reached
303             }
304             return data () [layout_type::element (k, size2_,
305                                                        l, lower_ + 1 + upper_)];
306 #else
307             // This is the netlib layout as described here: http://www.netlib.org/lapack/lug/node124.html
308             BOOST_UBLAS_CHECK( is_element_in_band( i, j) , bad_index());
309             return data () [hidden::banded_indexing<orientation_category>::get_index(size1_, size2_, lower_, upper_, i, j)];
310 #endif
311 
312         }
313 
314         // Element assignment
315         BOOST_UBLAS_INLINE
insert_element(size_type i,size_type j,const_reference t)316         reference insert_element (size_type i, size_type j, const_reference t) {
317             return (operator () (i, j) = t);
318         }
319         BOOST_UBLAS_INLINE
erase_element(size_type i,size_type j)320         void erase_element (size_type i, size_type j) {
321             operator () (i, j) = value_type/*zero*/();
322         }
323 
324         // Zeroing
325         BOOST_UBLAS_INLINE
clear()326         void clear () {
327             std::fill (data ().begin (), data ().end (), value_type/*zero*/());
328         }
329 
330         // Assignment
331         BOOST_UBLAS_INLINE
operator =(const banded_matrix & m)332         banded_matrix &operator = (const banded_matrix &m) {
333             size1_ = m.size1_;
334             size2_ = m.size2_;
335             lower_ = m.lower_;
336             upper_ = m.upper_;
337             data () = m.data ();
338             return *this;
339         }
340         BOOST_UBLAS_INLINE
assign_temporary(banded_matrix & m)341         banded_matrix &assign_temporary (banded_matrix &m) {
342             swap (m);
343             return *this;
344         }
345         template<class AE>
346         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)347         banded_matrix &operator = (const matrix_expression<AE> &ae) {
348             self_type temporary (ae, lower_, upper_);
349             return assign_temporary (temporary);
350         }
351         template<class AE>
352         BOOST_UBLAS_INLINE
assign(const matrix_expression<AE> & ae)353         banded_matrix &assign (const matrix_expression<AE> &ae) {
354             matrix_assign<scalar_assign> (*this, ae);
355             return *this;
356         }
357         template<class AE>
358         BOOST_UBLAS_INLINE
operator +=(const matrix_expression<AE> & ae)359         banded_matrix& operator += (const matrix_expression<AE> &ae) {
360             self_type temporary (*this + ae, lower_, upper_);
361             return assign_temporary (temporary);
362         }
363         template<class AE>
364         BOOST_UBLAS_INLINE
plus_assign(const matrix_expression<AE> & ae)365         banded_matrix &plus_assign (const matrix_expression<AE> &ae) {
366             matrix_assign<scalar_plus_assign> (*this, ae);
367             return *this;
368         }
369         template<class AE>
370         BOOST_UBLAS_INLINE
operator -=(const matrix_expression<AE> & ae)371         banded_matrix& operator -= (const matrix_expression<AE> &ae) {
372             self_type temporary (*this - ae, lower_, upper_);
373             return assign_temporary (temporary);
374         }
375         template<class AE>
376         BOOST_UBLAS_INLINE
minus_assign(const matrix_expression<AE> & ae)377         banded_matrix &minus_assign (const matrix_expression<AE> &ae) {
378             matrix_assign<scalar_minus_assign> (*this, ae);
379             return *this;
380         }
381         template<class AT>
382         BOOST_UBLAS_INLINE
operator *=(const AT & at)383         banded_matrix& operator *= (const AT &at) {
384             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
385             return *this;
386         }
387         template<class AT>
388         BOOST_UBLAS_INLINE
operator /=(const AT & at)389         banded_matrix& operator /= (const AT &at) {
390             matrix_assign_scalar<scalar_divides_assign> (*this, at);
391             return *this;
392         }
393 
394         // Swapping
395         BOOST_UBLAS_INLINE
swap(banded_matrix & m)396         void swap (banded_matrix &m) {
397             if (this != &m) {
398                 std::swap (size1_, m.size1_);
399                 std::swap (size2_, m.size2_);
400                 std::swap (lower_, m.lower_);
401                 std::swap (upper_, m.upper_);
402                 data ().swap (m.data ());
403             }
404         }
405         BOOST_UBLAS_INLINE
swap(banded_matrix & m1,banded_matrix & m2)406         friend void swap (banded_matrix &m1, banded_matrix &m2) {
407             m1.swap (m2);
408         }
409 
410         // Iterator types
411 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
412         typedef indexed_iterator1<self_type, packed_random_access_iterator_tag> iterator1;
413         typedef indexed_iterator2<self_type, packed_random_access_iterator_tag> iterator2;
414         typedef indexed_const_iterator1<self_type, packed_random_access_iterator_tag> const_iterator1;
415         typedef indexed_const_iterator2<self_type, packed_random_access_iterator_tag> const_iterator2;
416 #else
417         class const_iterator1;
418         class iterator1;
419         class const_iterator2;
420         class iterator2;
421 #endif
422         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
423         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
424         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
425         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
426 
427         // Element lookup
428         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j) const429         const_iterator1 find1 (int rank, size_type i, size_type j) const {
430             if (rank == 1) {
431                 size_type lower_i = (std::max) (difference_type (j - upper_), difference_type (0));
432                 i = (std::max) (i, lower_i);
433                 size_type upper_i = (std::min) (j + 1 + lower_, size1_);
434                 i = (std::min) (i, upper_i);
435             }
436             return const_iterator1 (*this, i, j);
437         }
438         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j)439         iterator1 find1 (int rank, size_type i, size_type j) {
440             if (rank == 1) {
441                 size_type lower_i = (std::max) (difference_type (j - upper_), difference_type (0));
442                 i = (std::max) (i, lower_i);
443                 size_type upper_i = (std::min) (j + 1 + lower_, size1_);
444                 i = (std::min) (i, upper_i);
445             }
446             return iterator1 (*this, i, j);
447         }
448         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j) const449         const_iterator2 find2 (int rank, size_type i, size_type j) const {
450             if (rank == 1) {
451                 size_type lower_j = (std::max) (difference_type (i - lower_), difference_type (0));
452                 j = (std::max) (j, lower_j);
453                 size_type upper_j = (std::min) (i + 1 + upper_, size2_);
454                 j = (std::min) (j, upper_j);
455             }
456             return const_iterator2 (*this, i, j);
457         }
458         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j)459         iterator2 find2 (int rank, size_type i, size_type j) {
460             if (rank == 1) {
461                 size_type lower_j = (std::max) (difference_type (i - lower_), difference_type (0));
462                 j = (std::max) (j, lower_j);
463                 size_type upper_j = (std::min) (i + 1 + upper_, size2_);
464                 j = (std::min) (j, upper_j);
465             }
466             return iterator2 (*this, i, j);
467         }
468 
469         // Iterators simply are indices.
470 
471 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
472         class const_iterator1:
473             public container_const_reference<banded_matrix>,
474             public random_access_iterator_base<packed_random_access_iterator_tag,
475                                                const_iterator1, value_type> {
476         public:
477             typedef typename banded_matrix::value_type value_type;
478             typedef typename banded_matrix::difference_type difference_type;
479             typedef typename banded_matrix::const_reference reference;
480             typedef const typename banded_matrix::pointer pointer;
481 
482             typedef const_iterator2 dual_iterator_type;
483             typedef const_reverse_iterator2 dual_reverse_iterator_type;
484 
485             // Construction and destruction
486             BOOST_UBLAS_INLINE
const_iterator1()487             const_iterator1 ():
488                 container_const_reference<self_type> (), it1_ (), it2_ () {}
489             BOOST_UBLAS_INLINE
const_iterator1(const self_type & m,size_type it1,size_type it2)490             const_iterator1 (const self_type &m, size_type it1, size_type it2):
491                 container_const_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
492             BOOST_UBLAS_INLINE
const_iterator1(const iterator1 & it)493             const_iterator1 (const iterator1 &it):
494                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
495 
496             // Arithmetic
497             BOOST_UBLAS_INLINE
operator ++()498             const_iterator1 &operator ++ () {
499                 ++ it1_;
500                 return *this;
501             }
502             BOOST_UBLAS_INLINE
operator --()503             const_iterator1 &operator -- () {
504                 -- it1_;
505                 return *this;
506             }
507             BOOST_UBLAS_INLINE
operator +=(difference_type n)508             const_iterator1 &operator += (difference_type n) {
509                 it1_ += n;
510                 return *this;
511             }
512             BOOST_UBLAS_INLINE
operator -=(difference_type n)513             const_iterator1 &operator -= (difference_type n) {
514                 it1_ -= n;
515                 return *this;
516             }
517             BOOST_UBLAS_INLINE
operator -(const const_iterator1 & it) const518             difference_type operator - (const const_iterator1 &it) const {
519                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
520                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
521                 return it1_ - it.it1_;
522             }
523 
524             // Dereference
525             BOOST_UBLAS_INLINE
operator *() const526             const_reference operator * () const {
527                 return (*this) () (it1_, it2_);
528             }
529             BOOST_UBLAS_INLINE
operator [](difference_type n) const530             const_reference operator [] (difference_type n) const {
531                 return *(*this + n);
532             }
533 
534 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
535             BOOST_UBLAS_INLINE
536 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
537             typename self_type::
538 #endif
begin() const539             const_iterator2 begin () const {
540                 return (*this) ().find2 (1, it1_, 0);
541             }
542             BOOST_UBLAS_INLINE
543 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
544             typename self_type::
545 #endif
cbegin() const546             const_iterator2 cbegin () const {
547                 return begin ();
548             }
549             BOOST_UBLAS_INLINE
550 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
551             typename self_type::
552 #endif
end() const553             const_iterator2 end () const {
554                 return (*this) ().find2 (1, it1_, (*this) ().size2 ());
555             }
556             BOOST_UBLAS_INLINE
557 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
558             typename self_type::
559 #endif
cend() const560             const_iterator2 cend () const {
561                 return end ();
562             }
563             BOOST_UBLAS_INLINE
564 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
565             typename self_type::
566 #endif
rbegin() const567             const_reverse_iterator2 rbegin () const {
568                 return const_reverse_iterator2 (end ());
569             }
570             BOOST_UBLAS_INLINE
571 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
572             typename self_type::
573 #endif
crbegin() const574             const_reverse_iterator2 crbegin () const {
575                 return rbegin ();
576             }
577             BOOST_UBLAS_INLINE
578 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
579             typename self_type::
580 #endif
rend() const581             const_reverse_iterator2 rend () const {
582                 return const_reverse_iterator2 (begin ());
583             }
584             BOOST_UBLAS_INLINE
585 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
586             typename self_type::
587 #endif
crend() const588             const_reverse_iterator2 crend () const {
589                 return rend ();
590             }
591 #endif
592 
593             // Indices
594             BOOST_UBLAS_INLINE
index1() const595             size_type index1 () const {
596                 return it1_;
597             }
598             BOOST_UBLAS_INLINE
index2() const599             size_type index2 () const {
600                 return it2_;
601             }
602 
603             // Assignment
604             BOOST_UBLAS_INLINE
operator =(const const_iterator1 & it)605             const_iterator1 &operator = (const const_iterator1 &it) {
606                 container_const_reference<self_type>::assign (&it ());
607                 it1_ = it.it1_;
608                 it2_ = it.it2_;
609                 return *this;
610             }
611 
612             // Comparison
613             BOOST_UBLAS_INLINE
operator ==(const const_iterator1 & it) const614             bool operator == (const const_iterator1 &it) const {
615                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
616                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
617                 return it1_ == it.it1_;
618             }
619             BOOST_UBLAS_INLINE
operator <(const const_iterator1 & it) const620             bool operator < (const const_iterator1 &it) const {
621                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
622                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
623                 return it1_ < it.it1_;
624             }
625 
626         private:
627             size_type it1_;
628             size_type it2_;
629         };
630 #endif
631 
632         BOOST_UBLAS_INLINE
begin1() const633         const_iterator1 begin1 () const {
634             return find1 (0, 0, 0);
635         }
636         BOOST_UBLAS_INLINE
cbegin1() const637         const_iterator1 cbegin1 () const {
638             return begin1 ();
639         }
640         BOOST_UBLAS_INLINE
end1() const641         const_iterator1 end1 () const {
642             return find1 (0, size1_, 0);
643         }
644         BOOST_UBLAS_INLINE
cend1() const645         const_iterator1 cend1 () const {
646             return end1 ();
647         }
648 
649 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
650         class iterator1:
651             public container_reference<banded_matrix>,
652             public random_access_iterator_base<packed_random_access_iterator_tag,
653                                                iterator1, value_type> {
654         public:
655             typedef typename banded_matrix::value_type value_type;
656             typedef typename banded_matrix::difference_type difference_type;
657             typedef typename banded_matrix::reference reference;
658             typedef typename banded_matrix::pointer pointer;
659 
660             typedef iterator2 dual_iterator_type;
661             typedef reverse_iterator2 dual_reverse_iterator_type;
662 
663             // Construction and destruction
664             BOOST_UBLAS_INLINE
iterator1()665             iterator1 ():
666                 container_reference<self_type> (), it1_ (), it2_ () {}
667             BOOST_UBLAS_INLINE
iterator1(self_type & m,size_type it1,size_type it2)668             iterator1 (self_type &m, size_type it1, size_type it2):
669                 container_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
670 
671             // Arithmetic
672             BOOST_UBLAS_INLINE
operator ++()673             iterator1 &operator ++ () {
674                 ++ it1_;
675                 return *this;
676             }
677             BOOST_UBLAS_INLINE
operator --()678             iterator1 &operator -- () {
679                 -- it1_;
680                 return *this;
681             }
682             BOOST_UBLAS_INLINE
operator +=(difference_type n)683             iterator1 &operator += (difference_type n) {
684                 it1_ += n;
685                 return *this;
686             }
687             BOOST_UBLAS_INLINE
operator -=(difference_type n)688             iterator1 &operator -= (difference_type n) {
689                 it1_ -= n;
690                 return *this;
691             }
692             BOOST_UBLAS_INLINE
operator -(const iterator1 & it) const693             difference_type operator - (const iterator1 &it) const {
694                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
695                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
696                 return it1_ - it.it1_;
697             }
698 
699             // Dereference
700             BOOST_UBLAS_INLINE
operator *() const701             reference operator * () const {
702                 return (*this) ().at_element (it1_, it2_);
703             }
704             BOOST_UBLAS_INLINE
operator [](difference_type n) const705             reference operator [] (difference_type n) const {
706                 return *(*this + n);
707             }
708 
709 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
710             BOOST_UBLAS_INLINE
711 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
712             typename self_type::
713 #endif
begin() const714             iterator2 begin () const {
715                 return (*this) ().find2 (1, it1_, 0);
716             }
717             BOOST_UBLAS_INLINE
718 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
719             typename self_type::
720 #endif
end() const721             iterator2 end () const {
722                 return (*this) ().find2 (1, it1_, (*this) ().size2 ());
723             }
724 
725             BOOST_UBLAS_INLINE
726 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
727             typename self_type::
728 #endif
rbegin() const729             reverse_iterator2 rbegin () const {
730                 return reverse_iterator2 (end ());
731             }
732             BOOST_UBLAS_INLINE
733 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
734             typename self_type::
735 #endif
rend() const736             reverse_iterator2 rend () const {
737                 return reverse_iterator2 (begin ());
738             }
739 #endif
740 
741             // Indices
742             BOOST_UBLAS_INLINE
index1() const743             size_type index1 () const {
744                 return it1_;
745             }
746             BOOST_UBLAS_INLINE
index2() const747             size_type index2 () const {
748                 return it2_;
749             }
750 
751             // Assignment
752             BOOST_UBLAS_INLINE
operator =(const iterator1 & it)753             iterator1 &operator = (const iterator1 &it) {
754                 container_reference<self_type>::assign (&it ());
755                 it1_ = it.it1_;
756                 it2_ = it.it2_;
757                 return *this;
758             }
759 
760             // Comparison
761             BOOST_UBLAS_INLINE
operator ==(const iterator1 & it) const762             bool operator == (const iterator1 &it) const {
763                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
764                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
765                 return it1_ == it.it1_;
766             }
767             BOOST_UBLAS_INLINE
operator <(const iterator1 & it) const768             bool operator < (const iterator1 &it) const {
769                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
770                 BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
771                 return it1_ < it.it1_;
772             }
773 
774         private:
775             size_type it1_;
776             size_type it2_;
777 
778             friend class const_iterator1;
779         };
780 #endif
781 
782         BOOST_UBLAS_INLINE
begin1()783         iterator1 begin1 () {
784             return find1 (0, 0, 0);
785         }
786         BOOST_UBLAS_INLINE
end1()787         iterator1 end1 () {
788             return find1 (0, size1_, 0);
789         }
790 
791 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
792         class const_iterator2:
793             public container_const_reference<banded_matrix>,
794             public random_access_iterator_base<packed_random_access_iterator_tag,
795                                                const_iterator2, value_type> {
796         public:
797             typedef typename banded_matrix::value_type value_type;
798             typedef typename banded_matrix::difference_type difference_type;
799             typedef typename banded_matrix::const_reference reference;
800             typedef const typename banded_matrix::pointer pointer;
801 
802             typedef const_iterator1 dual_iterator_type;
803             typedef const_reverse_iterator1 dual_reverse_iterator_type;
804 
805             // Construction and destruction
806             BOOST_UBLAS_INLINE
const_iterator2()807             const_iterator2 ():
808                 container_const_reference<self_type> (), it1_ (), it2_ () {}
809             BOOST_UBLAS_INLINE
const_iterator2(const self_type & m,size_type it1,size_type it2)810             const_iterator2 (const self_type &m, size_type it1, size_type it2):
811                 container_const_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
812             BOOST_UBLAS_INLINE
const_iterator2(const iterator2 & it)813             const_iterator2 (const iterator2 &it):
814                 container_const_reference<self_type> (it ()), it1_ (it.it1_), it2_ (it.it2_) {}
815 
816             // Arithmetic
817             BOOST_UBLAS_INLINE
operator ++()818             const_iterator2 &operator ++ () {
819                 ++ it2_;
820                 return *this;
821             }
822             BOOST_UBLAS_INLINE
operator --()823             const_iterator2 &operator -- () {
824                 -- it2_;
825                 return *this;
826             }
827             BOOST_UBLAS_INLINE
operator +=(difference_type n)828             const_iterator2 &operator += (difference_type n) {
829                 it2_ += n;
830                 return *this;
831             }
832             BOOST_UBLAS_INLINE
operator -=(difference_type n)833             const_iterator2 &operator -= (difference_type n) {
834                 it2_ -= n;
835                 return *this;
836             }
837             BOOST_UBLAS_INLINE
operator -(const const_iterator2 & it) const838             difference_type operator - (const const_iterator2 &it) const {
839                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
840                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
841                 return it2_ - it.it2_;
842             }
843 
844             // Dereference
845             BOOST_UBLAS_INLINE
operator *() const846             const_reference operator * () const {
847                 return (*this) () (it1_, it2_);
848             }
849             BOOST_UBLAS_INLINE
operator [](difference_type n) const850             const_reference operator [] (difference_type n) const {
851                 return *(*this + n);
852             }
853 
854 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
855             BOOST_UBLAS_INLINE
856 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
857             typename self_type::
858 #endif
begin() const859             const_iterator1 begin () const {
860                 return (*this) ().find1 (1, 0, it2_);
861             }
862             BOOST_UBLAS_INLINE
863 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
864             typename self_type::
865 #endif
cbegin() const866             const_iterator1 cbegin () const {
867                 return begin ();
868             }
869             BOOST_UBLAS_INLINE
870 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
871             typename self_type::
872 #endif
end() const873             const_iterator1 end () const {
874                 return (*this) ().find1 (1, (*this) ().size1 (), it2_);
875             }
876             BOOST_UBLAS_INLINE
877 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
878             typename self_type::
879 #endif
cend() const880             const_iterator1 cend () const {
881                 return end();
882             }
883 
884             BOOST_UBLAS_INLINE
885 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
886             typename self_type::
887 #endif
rbegin() const888             const_reverse_iterator1 rbegin () const {
889                 return const_reverse_iterator1 (end ());
890             }
891             BOOST_UBLAS_INLINE
892 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
893             typename self_type::
894 #endif
crbegin() const895             const_reverse_iterator1 crbegin () const {
896                 return rbegin ();
897             }
898             BOOST_UBLAS_INLINE
899 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
900             typename self_type::
901 #endif
rend() const902             const_reverse_iterator1 rend () const {
903                 return const_reverse_iterator1 (begin ());
904             }
905             BOOST_UBLAS_INLINE
906 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
907             typename self_type::
908 #endif
crend() const909             const_reverse_iterator1 crend () const {
910                 return rend ();
911             }
912 #endif
913 
914             // Indices
915             BOOST_UBLAS_INLINE
index1() const916             size_type index1 () const {
917                 return it1_;
918             }
919             BOOST_UBLAS_INLINE
index2() const920             size_type index2 () const {
921                 return it2_;
922             }
923 
924             // Assignment
925             BOOST_UBLAS_INLINE
operator =(const const_iterator2 & it)926             const_iterator2 &operator = (const const_iterator2 &it) {
927                 container_const_reference<self_type>::assign (&it ());
928                 it1_ = it.it1_;
929                 it2_ = it.it2_;
930                 return *this;
931             }
932 
933             // Comparison
934             BOOST_UBLAS_INLINE
operator ==(const const_iterator2 & it) const935             bool operator == (const const_iterator2 &it) const {
936                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
937                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
938                 return it2_ == it.it2_;
939             }
940             BOOST_UBLAS_INLINE
operator <(const const_iterator2 & it) const941             bool operator < (const const_iterator2 &it) const {
942                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
943                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
944                 return it2_ < it.it2_;
945             }
946 
947         private:
948             size_type it1_;
949             size_type it2_;
950         };
951 #endif
952 
953         BOOST_UBLAS_INLINE
begin2() const954         const_iterator2 begin2 () const {
955             return find2 (0, 0, 0);
956         }
957         BOOST_UBLAS_INLINE
cbegin2() const958         const_iterator2 cbegin2 () const {
959             return begin2 ();
960         }
961         BOOST_UBLAS_INLINE
end2() const962         const_iterator2 end2 () const {
963             return find2 (0, 0, size2_);
964         }
965         BOOST_UBLAS_INLINE
cend2() const966         const_iterator2 cend2 () const {
967             return end2 ();
968         }
969 
970 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
971         class iterator2:
972             public container_reference<banded_matrix>,
973             public random_access_iterator_base<packed_random_access_iterator_tag,
974                                                iterator2, value_type> {
975         public:
976             typedef typename banded_matrix::value_type value_type;
977             typedef typename banded_matrix::difference_type difference_type;
978             typedef typename banded_matrix::reference reference;
979             typedef typename banded_matrix::pointer pointer;
980 
981             typedef iterator1 dual_iterator_type;
982             typedef reverse_iterator1 dual_reverse_iterator_type;
983 
984             // Construction and destruction
985             BOOST_UBLAS_INLINE
iterator2()986             iterator2 ():
987                 container_reference<self_type> (), it1_ (), it2_ () {}
988             BOOST_UBLAS_INLINE
iterator2(self_type & m,size_type it1,size_type it2)989             iterator2 (self_type &m, size_type it1, size_type it2):
990                 container_reference<self_type> (m), it1_ (it1), it2_ (it2) {}
991 
992             // Arithmetic
993             BOOST_UBLAS_INLINE
operator ++()994             iterator2 &operator ++ () {
995                 ++ it2_;
996                 return *this;
997             }
998             BOOST_UBLAS_INLINE
operator --()999             iterator2 &operator -- () {
1000                 -- it2_;
1001                 return *this;
1002             }
1003             BOOST_UBLAS_INLINE
operator +=(difference_type n)1004             iterator2 &operator += (difference_type n) {
1005                 it2_ += n;
1006                 return *this;
1007             }
1008             BOOST_UBLAS_INLINE
operator -=(difference_type n)1009             iterator2 &operator -= (difference_type n) {
1010                 it2_ -= n;
1011                 return *this;
1012             }
1013             BOOST_UBLAS_INLINE
operator -(const iterator2 & it) const1014             difference_type operator - (const iterator2 &it) const {
1015                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1016                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1017                 return it2_ - it.it2_;
1018             }
1019 
1020             // Dereference
1021             BOOST_UBLAS_INLINE
operator *() const1022             reference operator * () const {
1023                 return (*this) ().at_element (it1_, it2_);
1024             }
1025             BOOST_UBLAS_INLINE
operator [](difference_type n) const1026             reference operator [] (difference_type n) const {
1027                 return *(*this + n);
1028             }
1029 
1030 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1031             BOOST_UBLAS_INLINE
1032 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1033             typename self_type::
1034 #endif
begin() const1035             iterator1 begin () const {
1036                 return (*this) ().find1 (1, 0, it2_);
1037             }
1038             BOOST_UBLAS_INLINE
1039 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1040             typename self_type::
1041 #endif
end() const1042             iterator1 end () const {
1043                 return (*this) ().find1 (1, (*this) ().size1 (), it2_);
1044             }
1045             BOOST_UBLAS_INLINE
1046 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1047             typename self_type::
1048 #endif
rbegin() const1049             reverse_iterator1 rbegin () const {
1050                 return reverse_iterator1 (end ());
1051             }
1052             BOOST_UBLAS_INLINE
1053 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1054             typename self_type::
1055 #endif
rend() const1056             reverse_iterator1 rend () const {
1057                 return reverse_iterator1 (begin ());
1058             }
1059 #endif
1060 
1061             // Indices
1062             BOOST_UBLAS_INLINE
index1() const1063             size_type index1 () const {
1064                 return it1_;
1065             }
1066             BOOST_UBLAS_INLINE
index2() const1067             size_type index2 () const {
1068                 return it2_;
1069             }
1070 
1071             // Assignment
1072             BOOST_UBLAS_INLINE
operator =(const iterator2 & it)1073             iterator2 &operator = (const iterator2 &it) {
1074                 container_reference<self_type>::assign (&it ());
1075                 it1_ = it.it1_;
1076                 it2_ = it.it2_;
1077                 return *this;
1078             }
1079 
1080             // Comparison
1081             BOOST_UBLAS_INLINE
operator ==(const iterator2 & it) const1082             bool operator == (const iterator2 &it) const {
1083                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1084                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1085                 return it2_ == it.it2_;
1086             }
1087             BOOST_UBLAS_INLINE
operator <(const iterator2 & it) const1088             bool operator < (const iterator2 &it) const {
1089                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1090                 BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1091                 return it2_ < it.it2_;
1092             }
1093 
1094         private:
1095             size_type it1_;
1096             size_type it2_;
1097 
1098             friend class const_iterator2;
1099         };
1100 #endif
1101 
1102         BOOST_UBLAS_INLINE
begin2()1103         iterator2 begin2 () {
1104             return find2 (0, 0, 0);
1105         }
1106         BOOST_UBLAS_INLINE
end2()1107         iterator2 end2 () {
1108             return find2 (0, 0, size2_);
1109         }
1110 
1111         // Reverse iterators
1112 
1113         BOOST_UBLAS_INLINE
rbegin1() const1114         const_reverse_iterator1 rbegin1 () const {
1115             return const_reverse_iterator1 (end1 ());
1116         }
1117         BOOST_UBLAS_INLINE
crbegin1() const1118         const_reverse_iterator1 crbegin1 () const {
1119             return rbegin1 ();
1120         }
1121         BOOST_UBLAS_INLINE
rend1() const1122         const_reverse_iterator1 rend1 () const {
1123             return const_reverse_iterator1 (begin1 ());
1124         }
1125         BOOST_UBLAS_INLINE
crend1() const1126         const_reverse_iterator1 crend1 () const {
1127             return rend1 ();
1128         }
1129 
1130         BOOST_UBLAS_INLINE
rbegin1()1131         reverse_iterator1 rbegin1 () {
1132             return reverse_iterator1 (end1 ());
1133         }
1134         BOOST_UBLAS_INLINE
rend1()1135         reverse_iterator1 rend1 () {
1136             return reverse_iterator1 (begin1 ());
1137         }
1138 
1139         BOOST_UBLAS_INLINE
rbegin2() const1140         const_reverse_iterator2 rbegin2 () const {
1141             return const_reverse_iterator2 (end2 ());
1142         }
1143         BOOST_UBLAS_INLINE
crbegin2() const1144         const_reverse_iterator2 crbegin2 () const {
1145             return rbegin2 ();
1146         }
1147         BOOST_UBLAS_INLINE
rend2() const1148         const_reverse_iterator2 rend2 () const {
1149             return const_reverse_iterator2 (begin2 ());
1150         }
1151         BOOST_UBLAS_INLINE
crend2() const1152         const_reverse_iterator2 crend2 () const {
1153             return rend2 ();
1154         }
1155 
1156         BOOST_UBLAS_INLINE
rbegin2()1157         reverse_iterator2 rbegin2 () {
1158             return reverse_iterator2 (end2 ());
1159         }
1160         BOOST_UBLAS_INLINE
rend2()1161         reverse_iterator2 rend2 () {
1162             return reverse_iterator2 (begin2 ());
1163         }
1164 
1165     private:
1166         size_type size1_;
1167         size_type size2_;
1168         size_type lower_;
1169         size_type upper_;
1170         array_type data_;
1171         typedef const value_type const_value_type;
1172         static const_value_type zero_;
1173     };
1174 
1175     template<class T, class L, class A>
1176     typename banded_matrix<T, L, A>::const_value_type banded_matrix<T, L, A>::zero_ = value_type/*zero*/();
1177 
1178 
1179     /** \brief A diagonal matrix of values of type \c T, which is a specialization of a banded matrix
1180      *
1181      * For a \f$(m\times m)\f$-dimensional diagonal matrix, \f$0 \leq i < m\f$ and \f$0 \leq j < m\f$,
1182      * if \f$i\neq j\f$ then \f$b_{i,j}=0\f$. The default storage for diagonal matrices is packed.
1183      * Orientation and storage can also be specified. Default is \c row major \c unbounded_array.
1184      *
1185      * As a specialization of a banded matrix, the constructor of the diagonal matrix creates
1186      * a banded matrix with 0 upper and lower diagonals around the main diagonal and the matrix is
1187      * obviously a square matrix. Operations are optimized based on these 2 assumptions. It is
1188      * \b not required by the storage to initialize elements of the matrix.
1189      *
1190      * \tparam T the type of object stored in the matrix (like double, float, complex, etc...)
1191      * \tparam L the storage organization. It can be either \c row_major or \c column_major. Default is \c row_major
1192      * \tparam A the type of Storage array. Default is \c unbounded_array
1193      */
1194     template<class T, class L, class A>
1195     class diagonal_matrix:
1196         public banded_matrix<T, L, A> {
1197     public:
1198         typedef typename A::size_type size_type;
1199         typedef banded_matrix<T, L, A> matrix_type;
1200         typedef A array_type;
1201 
1202         // Construction and destruction
1203         BOOST_UBLAS_INLINE
diagonal_matrix()1204         diagonal_matrix ():
1205             matrix_type () {}
1206         BOOST_UBLAS_INLINE
diagonal_matrix(size_type size)1207         diagonal_matrix (size_type size):
1208             matrix_type (size, size) {}
1209         BOOST_UBLAS_INLINE
diagonal_matrix(size_type size,const array_type & data)1210         diagonal_matrix (size_type size, const array_type& data):
1211             matrix_type (size, size, 0, 0, data) {}
1212         BOOST_UBLAS_INLINE
diagonal_matrix(size_type size1,size_type size2)1213         diagonal_matrix (size_type size1, size_type size2):
1214             matrix_type (size1, size2) {}
1215         template<class AE>
1216         BOOST_UBLAS_INLINE
diagonal_matrix(const matrix_expression<AE> & ae)1217         diagonal_matrix (const matrix_expression<AE> &ae):
1218             matrix_type (ae) {}
1219         BOOST_UBLAS_INLINE
~diagonal_matrix()1220         ~diagonal_matrix () {}
1221 
1222         // Assignment
1223         BOOST_UBLAS_INLINE
operator =(const diagonal_matrix & m)1224         diagonal_matrix &operator = (const diagonal_matrix &m) {
1225             matrix_type::operator = (m);
1226             return *this;
1227         }
1228         template<class AE>
1229         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)1230         diagonal_matrix &operator = (const matrix_expression<AE> &ae) {
1231             matrix_type::operator = (ae);
1232             return *this;
1233         }
1234     };
1235 
1236     /** \brief A banded matrix adaptator: convert a any matrix into a banded matrix expression
1237      *
1238      * For a \f$(m\times n)\f$-dimensional matrix, the \c banded_adaptor will provide a banded matrix
1239      * with \f$l\f$ lower and \f$u\f$ upper diagonals and \f$0 \leq i < m\f$ and \f$0 \leq j < n\f$,
1240      * if \f$i>j+l\f$ or \f$i<j-u\f$ then \f$b_{i,j}=0\f$.
1241      *
1242      * Storage and location are based on those of the underlying matrix. This is important because
1243      * a \c banded_adaptor does not copy the matrix data to a new place. Therefore, modifying values
1244      * in a \c banded_adaptor matrix will also modify the underlying matrix too.
1245      *
1246      * \tparam M the type of matrix used to generate a banded matrix
1247      */
1248     template<class M>
1249     class banded_adaptor:
1250         public matrix_expression<banded_adaptor<M> > {
1251 
1252         typedef banded_adaptor<M> self_type;
1253 
1254     public:
1255 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1256         using matrix_expression<self_type>::operator ();
1257 #endif
1258         typedef const M const_matrix_type;
1259         typedef M matrix_type;
1260         typedef typename M::size_type size_type;
1261         typedef typename M::difference_type difference_type;
1262         typedef typename M::value_type value_type;
1263         typedef typename M::const_reference const_reference;
1264         typedef typename boost::mpl::if_<boost::is_const<M>,
1265                                           typename M::const_reference,
1266                                           typename M::reference>::type reference;
1267         typedef typename boost::mpl::if_<boost::is_const<M>,
1268                                           typename M::const_closure_type,
1269                                           typename M::closure_type>::type matrix_closure_type;
1270         typedef const self_type const_closure_type;
1271         typedef self_type closure_type;
1272         // Replaced by _temporary_traits to avoid type requirements on M
1273         //typedef typename M::vector_temporary_type vector_temporary_type;
1274         //typedef typename M::matrix_temporary_type matrix_temporary_type;
1275         typedef typename storage_restrict_traits<typename M::storage_category,
1276                                                  packed_proxy_tag>::storage_category storage_category;
1277         typedef typename M::orientation_category orientation_category;
1278 
1279         // Construction and destruction
1280         BOOST_UBLAS_INLINE
banded_adaptor(matrix_type & data,size_type lower=0,size_type upper=0)1281         banded_adaptor (matrix_type &data, size_type lower = 0, size_type upper = 0):
1282             matrix_expression<self_type> (),
1283             data_ (data), lower_ (lower), upper_ (upper) {}
1284         BOOST_UBLAS_INLINE
banded_adaptor(const banded_adaptor & m)1285         banded_adaptor (const banded_adaptor &m):
1286             matrix_expression<self_type> (),
1287             data_ (m.data_), lower_ (m.lower_), upper_ (m.upper_) {}
1288 
1289         // Accessors
1290         BOOST_UBLAS_INLINE
size1() const1291         size_type size1 () const {
1292             return data_.size1 ();
1293         }
1294         BOOST_UBLAS_INLINE
size2() const1295         size_type size2 () const {
1296             return data_.size2 ();
1297         }
1298         BOOST_UBLAS_INLINE
lower() const1299         size_type lower () const {
1300             return lower_;
1301         }
1302         BOOST_UBLAS_INLINE
upper() const1303         size_type upper () const {
1304             return upper_;
1305         }
1306 
1307         // Storage accessors
1308         BOOST_UBLAS_INLINE
data() const1309         const matrix_closure_type &data () const {
1310             return data_;
1311         }
1312         BOOST_UBLAS_INLINE
data()1313         matrix_closure_type &data () {
1314             return data_;
1315         }
1316 
1317 #if !defined (BOOST_UBLAS_OWN_BANDED)||(BOOST_UBLAS_LEGACY_BANDED)
1318         BOOST_UBLAS_INLINE
is_element_in_band(size_type i,size_type j) const1319         bool is_element_in_band(size_type i, size_type j) const{
1320             //return (upper_+i >= j) && i <= std::min(size1() - 1, j + lower_); // We don't need to check if i is outside because it is checked anyway in the accessors.
1321             return (upper_+i >= j) && i <= ( j + lower_); // Essentially this band has "infinite" positive dimensions
1322         }
1323 #endif
1324 
1325         // Element access
1326 #ifndef BOOST_UBLAS_PROXY_CONST_MEMBER
1327         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const1328         const_reference operator () (size_type i, size_type j) const {
1329             BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1330             BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1331 #ifdef BOOST_UBLAS_OWN_BANDED
1332             size_type k = (std::max) (i, j);
1333             size_type l = lower_ + j - i;
1334             if (k < (std::max) (size1 (), size2 ()) &&
1335                 l < lower_ + 1 + upper_)
1336                 return data () (i, j);
1337 #elif BOOST_UBLAS_LEGACY_BANDED
1338             size_type k = j;
1339             size_type l = upper_ + i - j;
1340             if (k < size2 () &&
1341                 l < lower_ + 1 + upper_)
1342                 return data () (i, j);
1343 #else
1344             if (is_element_in_band( i, j))
1345                 return data () (i, j);
1346 #endif
1347             return zero_;
1348         }
1349         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j)1350         reference operator () (size_type i, size_type j) {
1351             BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1352             BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1353 #ifdef BOOST_UBLAS_OWN_BANDED
1354             size_type k = (std::max) (i, j);
1355             size_type l = lower_ + j - i;
1356             if (k < (std::max) (size1 (), size2 ()) &&
1357                 l < lower_ + 1 + upper_)
1358                 return data () (i, j);
1359 #elif BOOST_UBLAS_LEGACY_BANDED
1360             size_type k = j;
1361             size_type l = upper_ + i - j;
1362             if (k < size2 () &&
1363                 l < lower_ + 1 + upper_)
1364                 return data () (i, j);
1365 #else
1366             if (is_element_in_band( i, j))
1367                 return data () (i, j);
1368 #endif
1369 #ifndef BOOST_UBLAS_REFERENCE_CONST_MEMBER
1370             bad_index ().raise ();
1371 #endif
1372             return const_cast<reference>(zero_);
1373         }
1374 #else
1375         BOOST_UBLAS_INLINE
operator ()(size_type i,size_type j) const1376         reference operator () (size_type i, size_type j) const {
1377             BOOST_UBLAS_CHECK (i < size1 (), bad_index ());
1378             BOOST_UBLAS_CHECK (j < size2 (), bad_index ());
1379 #ifdef BOOST_UBLAS_OWN_BANDED
1380             size_type k = (std::max) (i, j);
1381             size_type l = lower_ + j - i;
1382             if (k < (std::max) (size1 (), size2 ()) &&
1383                 l < lower_ + 1 + upper_)
1384                 return data () (i, j);
1385 #elif BOOST_UBLAS_LEGACY_BANDED
1386             size_type k = j;
1387             size_type l = upper_ + i - j;
1388             if (k < size2 () &&
1389                 l < lower_ + 1 + upper_)
1390                 return data () (i, j);
1391 #else
1392             if (is_element_in_band( i, j))
1393                 return data () (i, j);
1394 #endif
1395 #ifndef BOOST_UBLAS_REFERENCE_CONST_MEMBER
1396             bad_index ().raise ();
1397 #endif
1398             return const_cast<reference>(zero_);
1399         }
1400 #endif
1401 
1402         // Assignment
1403         BOOST_UBLAS_INLINE
operator =(const banded_adaptor & m)1404         banded_adaptor &operator = (const banded_adaptor &m) {
1405             matrix_assign<scalar_assign> (*this, m);
1406             return *this;
1407         }
1408         BOOST_UBLAS_INLINE
assign_temporary(banded_adaptor & m)1409         banded_adaptor &assign_temporary (banded_adaptor &m) {
1410             *this = m;
1411             return *this;
1412         }
1413         template<class AE>
1414         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)1415         banded_adaptor &operator = (const matrix_expression<AE> &ae) {
1416             matrix_assign<scalar_assign> (*this, matrix<value_type> (ae));
1417             return *this;
1418         }
1419         template<class AE>
1420         BOOST_UBLAS_INLINE
assign(const matrix_expression<AE> & ae)1421         banded_adaptor &assign (const matrix_expression<AE> &ae) {
1422             matrix_assign<scalar_assign> (*this, ae);
1423             return *this;
1424         }
1425         template<class AE>
1426         BOOST_UBLAS_INLINE
operator +=(const matrix_expression<AE> & ae)1427         banded_adaptor& operator += (const matrix_expression<AE> &ae) {
1428             matrix_assign<scalar_assign> (*this, matrix<value_type> (*this + ae));
1429             return *this;
1430         }
1431         template<class AE>
1432         BOOST_UBLAS_INLINE
plus_assign(const matrix_expression<AE> & ae)1433         banded_adaptor &plus_assign (const matrix_expression<AE> &ae) {
1434             matrix_assign<scalar_plus_assign> (*this, ae);
1435             return *this;
1436         }
1437         template<class AE>
1438         BOOST_UBLAS_INLINE
operator -=(const matrix_expression<AE> & ae)1439         banded_adaptor& operator -= (const matrix_expression<AE> &ae) {
1440             matrix_assign<scalar_assign> (*this, matrix<value_type> (*this - ae));
1441             return *this;
1442         }
1443         template<class AE>
1444         BOOST_UBLAS_INLINE
minus_assign(const matrix_expression<AE> & ae)1445         banded_adaptor &minus_assign (const matrix_expression<AE> &ae) {
1446             matrix_assign<scalar_minus_assign> (*this, ae);
1447             return *this;
1448         }
1449         template<class AT>
1450         BOOST_UBLAS_INLINE
operator *=(const AT & at)1451         banded_adaptor& operator *= (const AT &at) {
1452             matrix_assign_scalar<scalar_multiplies_assign> (*this, at);
1453             return *this;
1454         }
1455         template<class AT>
1456         BOOST_UBLAS_INLINE
operator /=(const AT & at)1457         banded_adaptor& operator /= (const AT &at) {
1458             matrix_assign_scalar<scalar_divides_assign> (*this, at);
1459             return *this;
1460         }
1461 
1462         // Closure comparison
1463         BOOST_UBLAS_INLINE
same_closure(const banded_adaptor & ba) const1464         bool same_closure (const banded_adaptor &ba) const {
1465             return (*this).data ().same_closure (ba.data ());
1466         }
1467 
1468         // Swapping
1469         BOOST_UBLAS_INLINE
swap(banded_adaptor & m)1470         void swap (banded_adaptor &m) {
1471             if (this != &m) {
1472                 BOOST_UBLAS_CHECK (lower_ == m.lower_, bad_size ());
1473                 BOOST_UBLAS_CHECK (upper_ == m.upper_, bad_size ());
1474                 matrix_swap<scalar_swap> (*this, m);
1475             }
1476         }
1477         BOOST_UBLAS_INLINE
swap(banded_adaptor & m1,banded_adaptor & m2)1478         friend void swap (banded_adaptor &m1, banded_adaptor &m2) {
1479             m1.swap (m2);
1480         }
1481 
1482         // Iterator types
1483     private:
1484         // Use the matrix iterator
1485         typedef typename M::const_iterator1 const_subiterator1_type;
1486         typedef typename boost::mpl::if_<boost::is_const<M>,
1487                                           typename M::const_iterator1,
1488                                           typename M::iterator1>::type subiterator1_type;
1489         typedef typename M::const_iterator2 const_subiterator2_type;
1490         typedef typename boost::mpl::if_<boost::is_const<M>,
1491                                           typename M::const_iterator2,
1492                                           typename M::iterator2>::type subiterator2_type;
1493 
1494     public:
1495 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1496         typedef indexed_iterator1<self_type, packed_random_access_iterator_tag> iterator1;
1497         typedef indexed_iterator2<self_type, packed_random_access_iterator_tag> iterator2;
1498         typedef indexed_const_iterator1<self_type, packed_random_access_iterator_tag> const_iterator1;
1499         typedef indexed_const_iterator2<self_type, packed_random_access_iterator_tag> const_iterator2;
1500 #else
1501         class const_iterator1;
1502         class iterator1;
1503         class const_iterator2;
1504         class iterator2;
1505 #endif
1506         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
1507         typedef reverse_iterator_base1<iterator1> reverse_iterator1;
1508         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
1509         typedef reverse_iterator_base2<iterator2> reverse_iterator2;
1510 
1511         // Element lookup
1512         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j) const1513         const_iterator1 find1 (int rank, size_type i, size_type j) const {
1514             if (rank == 1) {
1515                 size_type lower_i = (std::max) (difference_type (j - upper_), difference_type (0));
1516                 i = (std::max) (i, lower_i);
1517                 size_type upper_i = (std::min) (j + 1 + lower_, size1 ());
1518                 i = (std::min) (i, upper_i);
1519             }
1520             return const_iterator1 (*this, data ().find1 (rank, i, j));
1521         }
1522         BOOST_UBLAS_INLINE
find1(int rank,size_type i,size_type j)1523         iterator1 find1 (int rank, size_type i, size_type j) {
1524             if (rank == 1) {
1525                 size_type lower_i = (std::max) (difference_type (j - upper_), difference_type (0));
1526                 i = (std::max) (i, lower_i);
1527                 size_type upper_i = (std::min) (j + 1 + lower_, size1 ());
1528                 i = (std::min) (i, upper_i);
1529             }
1530             return iterator1 (*this, data ().find1 (rank, i, j));
1531         }
1532         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j) const1533         const_iterator2 find2 (int rank, size_type i, size_type j) const {
1534             if (rank == 1) {
1535                 size_type lower_j = (std::max) (difference_type (i - lower_), difference_type (0));
1536                 j = (std::max) (j, lower_j);
1537                 size_type upper_j = (std::min) (i + 1 + upper_, size2 ());
1538                 j = (std::min) (j, upper_j);
1539             }
1540             return const_iterator2 (*this, data ().find2 (rank, i, j));
1541         }
1542         BOOST_UBLAS_INLINE
find2(int rank,size_type i,size_type j)1543         iterator2 find2 (int rank, size_type i, size_type j) {
1544             if (rank == 1) {
1545                 size_type lower_j = (std::max) (difference_type (i - lower_), difference_type (0));
1546                 j = (std::max) (j, lower_j);
1547                 size_type upper_j = (std::min) (i + 1 + upper_, size2 ());
1548                 j = (std::min) (j, upper_j);
1549             }
1550             return iterator2 (*this, data ().find2 (rank, i, j));
1551         }
1552 
1553         // Iterators simply are indices.
1554 
1555 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1556         class const_iterator1:
1557             public container_const_reference<banded_adaptor>,
1558             public random_access_iterator_base<typename iterator_restrict_traits<
1559                                                    typename const_subiterator1_type::iterator_category, packed_random_access_iterator_tag>::iterator_category,
1560                                                const_iterator1, value_type> {
1561         public:
1562             typedef typename const_subiterator1_type::value_type value_type;
1563             typedef typename const_subiterator1_type::difference_type difference_type;
1564             typedef typename const_subiterator1_type::reference reference;
1565             typedef typename const_subiterator1_type::pointer pointer;
1566 
1567             typedef const_iterator2 dual_iterator_type;
1568             typedef const_reverse_iterator2 dual_reverse_iterator_type;
1569 
1570             // Construction and destruction
1571             BOOST_UBLAS_INLINE
const_iterator1()1572             const_iterator1 ():
1573                 container_const_reference<self_type> (), it1_ () {}
1574             BOOST_UBLAS_INLINE
const_iterator1(const self_type & m,const const_subiterator1_type & it1)1575             const_iterator1 (const self_type &m, const const_subiterator1_type &it1):
1576                 container_const_reference<self_type> (m), it1_ (it1) {}
1577             BOOST_UBLAS_INLINE
const_iterator1(const iterator1 & it)1578             const_iterator1 (const iterator1 &it):
1579                 container_const_reference<self_type> (it ()), it1_ (it.it1_) {}
1580 
1581             // Arithmetic
1582             BOOST_UBLAS_INLINE
operator ++()1583             const_iterator1 &operator ++ () {
1584                 ++ it1_;
1585                 return *this;
1586             }
1587             BOOST_UBLAS_INLINE
operator --()1588             const_iterator1 &operator -- () {
1589                 -- it1_;
1590                 return *this;
1591             }
1592             BOOST_UBLAS_INLINE
operator +=(difference_type n)1593             const_iterator1 &operator += (difference_type n) {
1594                 it1_ += n;
1595                 return *this;
1596             }
1597             BOOST_UBLAS_INLINE
operator -=(difference_type n)1598             const_iterator1 &operator -= (difference_type n) {
1599                 it1_ -= n;
1600                 return *this;
1601             }
1602             BOOST_UBLAS_INLINE
operator -(const const_iterator1 & it) const1603             difference_type operator - (const const_iterator1 &it) const {
1604                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1605                 return it1_ - it.it1_;
1606             }
1607 
1608             // Dereference
1609             BOOST_UBLAS_INLINE
operator *() const1610             const_reference operator * () const {
1611                 size_type i = index1 ();
1612                 size_type j = index2 ();
1613                 BOOST_UBLAS_CHECK (i < (*this) ().size1 (), bad_index ());
1614                 BOOST_UBLAS_CHECK (j < (*this) ().size2 (), bad_index ());
1615 #ifdef BOOST_UBLAS_OWN_BANDED
1616                 size_type k = (std::max) (i, j);
1617                 size_type l = (*this) ().lower () + j - i;
1618                 if (k < (std::max) ((*this) ().size1 (), (*this) ().size2 ()) &&
1619                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1620                     return *it1_;
1621 #else
1622                 size_type k = j;
1623                 size_type l = (*this) ().upper () + i - j;
1624                 if (k < (*this) ().size2 () &&
1625                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1626                     return *it1_;
1627 #endif
1628                 return (*this) () (i, j);
1629             }
1630             BOOST_UBLAS_INLINE
operator [](difference_type n) const1631             const_reference operator [] (difference_type n) const {
1632                 return *(*this + n);
1633             }
1634 
1635 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1636             BOOST_UBLAS_INLINE
1637 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1638             typename self_type::
1639 #endif
begin() const1640             const_iterator2 begin () const {
1641                 return (*this) ().find2 (1, index1 (), 0);
1642             }
1643             BOOST_UBLAS_INLINE
1644 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1645             typename self_type::
1646 #endif
cbegin() const1647             const_iterator2 cbegin () const {
1648                 return begin ();
1649             }
1650             BOOST_UBLAS_INLINE
1651 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1652             typename self_type::
1653 #endif
end() const1654             const_iterator2 end () const {
1655                 return (*this) ().find2 (1, index1 (), (*this) ().size2 ());
1656             }
1657             BOOST_UBLAS_INLINE
1658 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1659             typename self_type::
1660 #endif
cend() const1661             const_iterator2 cend () const {
1662                 return end ();
1663             }
1664             BOOST_UBLAS_INLINE
1665 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1666             typename self_type::
1667 #endif
rbegin() const1668             const_reverse_iterator2 rbegin () const {
1669                 return const_reverse_iterator2 (end ());
1670             }
1671             BOOST_UBLAS_INLINE
1672 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1673             typename self_type::
1674 #endif
crbegin() const1675             const_reverse_iterator2 crbegin () const {
1676                 return rbegin ();
1677             }
1678             BOOST_UBLAS_INLINE
1679 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1680             typename self_type::
1681 #endif
rend() const1682             const_reverse_iterator2 rend () const {
1683                 return const_reverse_iterator2 (begin ());
1684             }
1685             BOOST_UBLAS_INLINE
1686 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1687             typename self_type::
1688 #endif
crend() const1689             const_reverse_iterator2 crend () const {
1690                 return rend ();
1691             }
1692 #endif
1693 
1694             // Indices
1695             BOOST_UBLAS_INLINE
index1() const1696             size_type index1 () const {
1697                 return it1_.index1 ();
1698             }
1699             BOOST_UBLAS_INLINE
index2() const1700             size_type index2 () const {
1701                 return it1_.index2 ();
1702             }
1703 
1704             // Assignment
1705             BOOST_UBLAS_INLINE
operator =(const const_iterator1 & it)1706             const_iterator1 &operator = (const const_iterator1 &it) {
1707                 container_const_reference<self_type>::assign (&it ());
1708                 it1_ = it.it1_;
1709                 return *this;
1710             }
1711 
1712             // Comparison
1713             BOOST_UBLAS_INLINE
operator ==(const const_iterator1 & it) const1714             bool operator == (const const_iterator1 &it) const {
1715                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1716                 return it1_ == it.it1_;
1717             }
1718             BOOST_UBLAS_INLINE
operator <(const const_iterator1 & it) const1719             bool operator < (const const_iterator1 &it) const {
1720                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1721                 return it1_ < it.it1_;
1722             }
1723 
1724         private:
1725             const_subiterator1_type it1_;
1726         };
1727 #endif
1728 
1729         BOOST_UBLAS_INLINE
begin1() const1730         const_iterator1 begin1 () const {
1731             return find1 (0, 0, 0);
1732         }
1733         BOOST_UBLAS_INLINE
cbegin1() const1734         const_iterator1 cbegin1 () const {
1735             return begin1 ();
1736         }
1737         BOOST_UBLAS_INLINE
end1() const1738         const_iterator1 end1 () const {
1739             return find1 (0, size1 (), 0);
1740         }
1741         BOOST_UBLAS_INLINE
cend1() const1742         const_iterator1 cend1 () const {
1743             return end1 ();
1744         }
1745 
1746 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1747         class iterator1:
1748             public container_reference<banded_adaptor>,
1749             public random_access_iterator_base<typename iterator_restrict_traits<
1750                                                    typename subiterator1_type::iterator_category, packed_random_access_iterator_tag>::iterator_category,
1751                                                iterator1, value_type> {
1752         public:
1753             typedef typename subiterator1_type::value_type value_type;
1754             typedef typename subiterator1_type::difference_type difference_type;
1755             typedef typename subiterator1_type::reference reference;
1756             typedef typename subiterator1_type::pointer pointer;
1757 
1758             typedef iterator2 dual_iterator_type;
1759             typedef reverse_iterator2 dual_reverse_iterator_type;
1760 
1761             // Construction and destruction
1762             BOOST_UBLAS_INLINE
iterator1()1763             iterator1 ():
1764                 container_reference<self_type> (), it1_ () {}
1765             BOOST_UBLAS_INLINE
iterator1(self_type & m,const subiterator1_type & it1)1766             iterator1 (self_type &m, const subiterator1_type &it1):
1767                 container_reference<self_type> (m), it1_ (it1) {}
1768 
1769             // Arithmetic
1770             BOOST_UBLAS_INLINE
operator ++()1771             iterator1 &operator ++ () {
1772                 ++ it1_;
1773                 return *this;
1774             }
1775             BOOST_UBLAS_INLINE
operator --()1776             iterator1 &operator -- () {
1777                 -- it1_;
1778                 return *this;
1779             }
1780             BOOST_UBLAS_INLINE
operator +=(difference_type n)1781             iterator1 &operator += (difference_type n) {
1782                 it1_ += n;
1783                 return *this;
1784             }
1785             BOOST_UBLAS_INLINE
operator -=(difference_type n)1786             iterator1 &operator -= (difference_type n) {
1787                 it1_ -= n;
1788                 return *this;
1789             }
1790             BOOST_UBLAS_INLINE
operator -(const iterator1 & it) const1791             difference_type operator - (const iterator1 &it) const {
1792                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1793                 return it1_ - it.it1_;
1794             }
1795 
1796             // Dereference
1797             BOOST_UBLAS_INLINE
operator *() const1798             reference operator * () const {
1799                 size_type i = index1 ();
1800                 size_type j = index2 ();
1801                 BOOST_UBLAS_CHECK (i < (*this) ().size1 (), bad_index ());
1802                 BOOST_UBLAS_CHECK (j < (*this) ().size2 (), bad_index ());
1803 #ifdef BOOST_UBLAS_OWN_BANDED
1804                 size_type k = (std::max) (i, j);
1805                 size_type l = (*this) ().lower () + j - i;
1806                 if (k < (std::max) ((*this) ().size1 (), (*this) ().size2 ()) &&
1807                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1808                     return *it1_;
1809 #else
1810                 size_type k = j;
1811                 size_type l = (*this) ().upper () + i - j;
1812                 if (k < (*this) ().size2 () &&
1813                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1814                     return *it1_;
1815 #endif
1816                 return (*this) () (i, j);
1817             }
1818             BOOST_UBLAS_INLINE
operator [](difference_type n) const1819             reference operator [] (difference_type n) const {
1820                 return *(*this + n);
1821             }
1822 
1823 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1824             BOOST_UBLAS_INLINE
1825 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1826             typename self_type::
1827 #endif
begin() const1828             iterator2 begin () const {
1829                 return (*this) ().find2 (1, index1 (), 0);
1830             }
1831             BOOST_UBLAS_INLINE
1832 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1833             typename self_type::
1834 #endif
end() const1835             iterator2 end () const {
1836                 return (*this) ().find2 (1, index1 (), (*this) ().size2 ());
1837             }
1838             BOOST_UBLAS_INLINE
1839 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1840             typename self_type::
1841 #endif
rbegin() const1842             reverse_iterator2 rbegin () const {
1843                 return reverse_iterator2 (end ());
1844             }
1845             BOOST_UBLAS_INLINE
1846 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1847             typename self_type::
1848 #endif
rend() const1849             reverse_iterator2 rend () const {
1850                 return reverse_iterator2 (begin ());
1851             }
1852 #endif
1853 
1854             // Indices
1855             BOOST_UBLAS_INLINE
index1() const1856             size_type index1 () const {
1857                 return it1_.index1 ();
1858             }
1859             BOOST_UBLAS_INLINE
index2() const1860             size_type index2 () const {
1861                 return it1_.index2 ();
1862             }
1863 
1864             // Assignment
1865             BOOST_UBLAS_INLINE
operator =(const iterator1 & it)1866             iterator1 &operator = (const iterator1 &it) {
1867                 container_reference<self_type>::assign (&it ());
1868                 it1_ = it.it1_;
1869                 return *this;
1870             }
1871 
1872             // Comparison
1873             BOOST_UBLAS_INLINE
operator ==(const iterator1 & it) const1874             bool operator == (const iterator1 &it) const {
1875                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1876                 return it1_ == it.it1_;
1877             }
1878             BOOST_UBLAS_INLINE
operator <(const iterator1 & it) const1879             bool operator < (const iterator1 &it) const {
1880                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1881                 return it1_ < it.it1_;
1882             }
1883 
1884         private:
1885             subiterator1_type it1_;
1886 
1887             friend class const_iterator1;
1888         };
1889 #endif
1890 
1891         BOOST_UBLAS_INLINE
begin1()1892         iterator1 begin1 () {
1893             return find1 (0, 0, 0);
1894         }
1895         BOOST_UBLAS_INLINE
end1()1896         iterator1 end1 () {
1897             return find1 (0, size1 (), 0);
1898         }
1899 
1900 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1901         class const_iterator2:
1902             public container_const_reference<banded_adaptor>,
1903             public random_access_iterator_base<packed_random_access_iterator_tag,
1904                                                const_iterator2, value_type> {
1905         public:
1906             typedef typename iterator_restrict_traits<typename const_subiterator2_type::iterator_category,
1907                                                       packed_random_access_iterator_tag>::iterator_category iterator_category;
1908             typedef typename const_subiterator2_type::value_type value_type;
1909             typedef typename const_subiterator2_type::difference_type difference_type;
1910             typedef typename const_subiterator2_type::reference reference;
1911             typedef typename const_subiterator2_type::pointer pointer;
1912 
1913             typedef const_iterator1 dual_iterator_type;
1914             typedef const_reverse_iterator1 dual_reverse_iterator_type;
1915 
1916             // Construction and destruction
1917             BOOST_UBLAS_INLINE
const_iterator2()1918             const_iterator2 ():
1919                 container_const_reference<self_type> (), it2_ () {}
1920             BOOST_UBLAS_INLINE
const_iterator2(const self_type & m,const const_subiterator2_type & it2)1921             const_iterator2 (const self_type &m, const const_subiterator2_type &it2):
1922                 container_const_reference<self_type> (m), it2_ (it2) {}
1923             BOOST_UBLAS_INLINE
const_iterator2(const iterator2 & it)1924             const_iterator2 (const iterator2 &it):
1925                 container_const_reference<self_type> (it ()), it2_ (it.it2_) {}
1926 
1927             // Arithmetic
1928             BOOST_UBLAS_INLINE
operator ++()1929             const_iterator2 &operator ++ () {
1930                 ++ it2_;
1931                 return *this;
1932             }
1933             BOOST_UBLAS_INLINE
operator --()1934             const_iterator2 &operator -- () {
1935                 -- it2_;
1936                 return *this;
1937             }
1938             BOOST_UBLAS_INLINE
operator +=(difference_type n)1939             const_iterator2 &operator += (difference_type n) {
1940                 it2_ += n;
1941                 return *this;
1942             }
1943             BOOST_UBLAS_INLINE
operator -=(difference_type n)1944             const_iterator2 &operator -= (difference_type n) {
1945                 it2_ -= n;
1946                 return *this;
1947             }
1948             BOOST_UBLAS_INLINE
operator -(const const_iterator2 & it) const1949             difference_type operator - (const const_iterator2 &it) const {
1950                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
1951                 return it2_ - it.it2_;
1952             }
1953 
1954             // Dereference
1955             BOOST_UBLAS_INLINE
operator *() const1956             const_reference operator * () const {
1957                 size_type i = index1 ();
1958                 size_type j = index2 ();
1959                 BOOST_UBLAS_CHECK (i < (*this) ().size1 (), bad_index ());
1960                 BOOST_UBLAS_CHECK (j < (*this) ().size2 (), bad_index ());
1961 #ifdef BOOST_UBLAS_OWN_BANDED
1962                 size_type k = (std::max) (i, j);
1963                 size_type l = (*this) ().lower () + j - i;
1964                 if (k < (std::max) ((*this) ().size1 (), (*this) ().size2 ()) &&
1965                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1966                     return *it2_;
1967 #else
1968                 size_type k = j;
1969                 size_type l = (*this) ().upper () + i - j;
1970                 if (k < (*this) ().size2 () &&
1971                     l < (*this) ().lower () + 1 + (*this) ().upper ())
1972                     return *it2_;
1973 #endif
1974                 return (*this) () (i, j);
1975             }
1976             BOOST_UBLAS_INLINE
operator [](difference_type n) const1977             const_reference operator [] (difference_type n) const {
1978                 return *(*this + n);
1979             }
1980 
1981 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1982             BOOST_UBLAS_INLINE
1983 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1984             typename self_type::
1985 #endif
begin() const1986             const_iterator1 begin () const {
1987                 return (*this) ().find1 (1, 0, index2 ());
1988             }
1989             BOOST_UBLAS_INLINE
1990 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1991             typename self_type::
1992 #endif
cbegin() const1993             const_iterator1 cbegin () const {
1994                 return begin ();
1995             }
1996             BOOST_UBLAS_INLINE
1997 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
1998             typename self_type::
1999 #endif
end() const2000             const_iterator1 end () const {
2001                 return (*this) ().find1 (1, (*this) ().size1 (), index2 ());
2002             }
2003             BOOST_UBLAS_INLINE
2004 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2005             typename self_type::
2006 #endif
cend() const2007             const_iterator1 cend () const {
2008                 return end ();
2009             }
2010             BOOST_UBLAS_INLINE
2011 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2012             typename self_type::
2013 #endif
rbegin() const2014             const_reverse_iterator1 rbegin () const {
2015                 return const_reverse_iterator1 (end ());
2016             }
2017             BOOST_UBLAS_INLINE
2018 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2019             typename self_type::
2020 #endif
crbegin() const2021             const_reverse_iterator1 crbegin () const {
2022                 return rbegin ();
2023             }
2024             BOOST_UBLAS_INLINE
2025 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2026             typename self_type::
2027 #endif
rend() const2028             const_reverse_iterator1 rend () const {
2029                 return const_reverse_iterator1 (begin ());
2030             }
2031             BOOST_UBLAS_INLINE
2032 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2033             typename self_type::
2034 #endif
crend() const2035             const_reverse_iterator1 crend () const {
2036                 return rend ();
2037             }
2038 #endif
2039 
2040             // Indices
2041             BOOST_UBLAS_INLINE
index1() const2042             size_type index1 () const {
2043                 return it2_.index1 ();
2044             }
2045             BOOST_UBLAS_INLINE
index2() const2046             size_type index2 () const {
2047                 return it2_.index2 ();
2048             }
2049 
2050             // Assignment
2051             BOOST_UBLAS_INLINE
operator =(const const_iterator2 & it)2052             const_iterator2 &operator = (const const_iterator2 &it) {
2053                 container_const_reference<self_type>::assign (&it ());
2054                 it2_ = it.it2_;
2055                 return *this;
2056             }
2057 
2058             // Comparison
2059             BOOST_UBLAS_INLINE
operator ==(const const_iterator2 & it) const2060             bool operator == (const const_iterator2 &it) const {
2061                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2062                 return it2_ == it.it2_;
2063             }
2064             BOOST_UBLAS_INLINE
operator <(const const_iterator2 & it) const2065             bool operator < (const const_iterator2 &it) const {
2066                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2067                 return it2_ < it.it2_;
2068             }
2069 
2070         private:
2071             const_subiterator2_type it2_;
2072         };
2073 #endif
2074 
2075         BOOST_UBLAS_INLINE
begin2() const2076         const_iterator2 begin2 () const {
2077             return find2 (0, 0, 0);
2078         }
2079         BOOST_UBLAS_INLINE
cbegin2() const2080         const_iterator2 cbegin2 () const {
2081             return begin2 ();
2082         }
2083         BOOST_UBLAS_INLINE
end2() const2084         const_iterator2 end2 () const {
2085             return find2 (0, 0, size2 ());
2086         }
2087         BOOST_UBLAS_INLINE
cend2() const2088         const_iterator2 cend2 () const {
2089             return end2 ();
2090         }
2091 
2092 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
2093         class iterator2:
2094             public container_reference<banded_adaptor>,
2095             public random_access_iterator_base<typename iterator_restrict_traits<
2096                                                    typename subiterator2_type::iterator_category, packed_random_access_iterator_tag>::iterator_category,
2097                                                iterator2, value_type> {
2098         public:
2099             typedef typename subiterator2_type::value_type value_type;
2100             typedef typename subiterator2_type::difference_type difference_type;
2101             typedef typename subiterator2_type::reference reference;
2102             typedef typename subiterator2_type::pointer pointer;
2103 
2104             typedef iterator1 dual_iterator_type;
2105             typedef reverse_iterator1 dual_reverse_iterator_type;
2106 
2107             // Construction and destruction
2108             BOOST_UBLAS_INLINE
iterator2()2109             iterator2 ():
2110                 container_reference<self_type> (), it2_ () {}
2111             BOOST_UBLAS_INLINE
iterator2(self_type & m,const subiterator2_type & it2)2112             iterator2 (self_type &m, const subiterator2_type &it2):
2113                 container_reference<self_type> (m), it2_ (it2) {}
2114 
2115             // Arithmetic
2116             BOOST_UBLAS_INLINE
operator ++()2117             iterator2 &operator ++ () {
2118                 ++ it2_;
2119                 return *this;
2120             }
2121             BOOST_UBLAS_INLINE
operator --()2122             iterator2 &operator -- () {
2123                 -- it2_;
2124                 return *this;
2125             }
2126             BOOST_UBLAS_INLINE
operator +=(difference_type n)2127             iterator2 &operator += (difference_type n) {
2128                 it2_ += n;
2129                 return *this;
2130             }
2131             BOOST_UBLAS_INLINE
operator -=(difference_type n)2132             iterator2 &operator -= (difference_type n) {
2133                 it2_ -= n;
2134                 return *this;
2135             }
2136             BOOST_UBLAS_INLINE
operator -(const iterator2 & it) const2137             difference_type operator - (const iterator2 &it) const {
2138                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2139                 return it2_ - it.it2_;
2140             }
2141 
2142             // Dereference
2143             BOOST_UBLAS_INLINE
operator *() const2144             reference operator * () const {
2145                 size_type i = index1 ();
2146                 size_type j = index2 ();
2147                 BOOST_UBLAS_CHECK (i < (*this) ().size1 (), bad_index ());
2148                 BOOST_UBLAS_CHECK (j < (*this) ().size2 (), bad_index ());
2149 #ifdef BOOST_UBLAS_OWN_BANDED
2150                 size_type k = (std::max) (i, j);
2151                 size_type l = (*this) ().lower () + j - i;
2152                 if (k < (std::max) ((*this) ().size1 (), (*this) ().size2 ()) &&
2153                     l < (*this) ().lower () + 1 + (*this) ().upper ())
2154                     return *it2_;
2155 #else
2156                 size_type k = j;
2157                 size_type l = (*this) ().upper () + i - j;
2158                 if (k < (*this) ().size2 () &&
2159                     l < (*this) ().lower () + 1 + (*this) ().upper ())
2160                     return *it2_;
2161 #endif
2162                 return (*this) () (i, j);
2163             }
2164             BOOST_UBLAS_INLINE
operator [](difference_type n) const2165             reference operator [] (difference_type n) const {
2166                 return *(*this + n);
2167             }
2168 
2169 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
2170             BOOST_UBLAS_INLINE
2171 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2172             typename self_type::
2173 #endif
begin() const2174             iterator1 begin () const {
2175                 return (*this) ().find1 (1, 0, index2 ());
2176             }
2177             BOOST_UBLAS_INLINE
2178 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2179             typename self_type::
2180 #endif
end() const2181             iterator1 end () const {
2182                 return (*this) ().find1 (1, (*this) ().size1 (), index2 ());
2183             }
2184             BOOST_UBLAS_INLINE
2185 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2186             typename self_type::
2187 #endif
rbegin() const2188             reverse_iterator1 rbegin () const {
2189                 return reverse_iterator1 (end ());
2190             }
2191             BOOST_UBLAS_INLINE
2192 #ifdef BOOST_UBLAS_MSVC_NESTED_CLASS_RELATION
2193             typename self_type::
2194 #endif
rend() const2195             reverse_iterator1 rend () const {
2196                 return reverse_iterator1 (begin ());
2197             }
2198 #endif
2199 
2200             // Indices
2201             BOOST_UBLAS_INLINE
index1() const2202             size_type index1 () const {
2203                 return it2_.index1 ();
2204             }
2205             BOOST_UBLAS_INLINE
index2() const2206             size_type index2 () const {
2207                 return it2_.index2 ();
2208             }
2209 
2210             // Assignment
2211             BOOST_UBLAS_INLINE
operator =(const iterator2 & it)2212             iterator2 &operator = (const iterator2 &it) {
2213                 container_reference<self_type>::assign (&it ());
2214                 it2_ = it.it2_;
2215                 return *this;
2216             }
2217 
2218             // Comparison
2219             BOOST_UBLAS_INLINE
operator ==(const iterator2 & it) const2220             bool operator == (const iterator2 &it) const {
2221                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2222                 return it2_ == it.it2_;
2223             }
2224             BOOST_UBLAS_INLINE
operator <(const iterator2 & it) const2225             bool operator < (const iterator2 &it) const {
2226                 BOOST_UBLAS_CHECK (&(*this) () == &it (), external_logic ());
2227                 return it2_ < it.it2_;
2228             }
2229 
2230         private:
2231             subiterator2_type it2_;
2232 
2233             friend class const_iterator2;
2234         };
2235 #endif
2236 
2237         BOOST_UBLAS_INLINE
begin2()2238         iterator2 begin2 () {
2239             return find2 (0, 0, 0);
2240         }
2241         BOOST_UBLAS_INLINE
end2()2242         iterator2 end2 () {
2243             return find2 (0, 0, size2 ());
2244         }
2245 
2246         // Reverse iterators
2247 
2248         BOOST_UBLAS_INLINE
rbegin1() const2249         const_reverse_iterator1 rbegin1 () const {
2250             return const_reverse_iterator1 (end1 ());
2251         }
2252         BOOST_UBLAS_INLINE
crbegin1() const2253         const_reverse_iterator1 crbegin1 () const {
2254             return rbegin1 ();
2255         }
2256         BOOST_UBLAS_INLINE
rend1() const2257         const_reverse_iterator1 rend1 () const {
2258             return const_reverse_iterator1 (begin1 ());
2259         }
2260         BOOST_UBLAS_INLINE
crend1() const2261         const_reverse_iterator1 crend1 () const {
2262             return rend1 ();
2263         }
2264 
2265         BOOST_UBLAS_INLINE
rbegin1()2266         reverse_iterator1 rbegin1 () {
2267             return reverse_iterator1 (end1 ());
2268         }
2269         BOOST_UBLAS_INLINE
rend1()2270         reverse_iterator1 rend1 () {
2271             return reverse_iterator1 (begin1 ());
2272         }
2273 
2274         BOOST_UBLAS_INLINE
rbegin2() const2275         const_reverse_iterator2 rbegin2 () const {
2276             return const_reverse_iterator2 (end2 ());
2277         }
2278         BOOST_UBLAS_INLINE
crbegin2() const2279         const_reverse_iterator2 crbegin2 () const {
2280             return rbegin2 ();
2281         }
2282         BOOST_UBLAS_INLINE
rend2() const2283         const_reverse_iterator2 rend2 () const {
2284             return const_reverse_iterator2 (begin2 ());
2285         }
2286         BOOST_UBLAS_INLINE
crend2() const2287         const_reverse_iterator2 crend2 () const {
2288             return rend2 ();
2289         }
2290 
2291         BOOST_UBLAS_INLINE
rbegin2()2292         reverse_iterator2 rbegin2 () {
2293             return reverse_iterator2 (end2 ());
2294         }
2295         BOOST_UBLAS_INLINE
rend2()2296         reverse_iterator2 rend2 () {
2297             return reverse_iterator2 (begin2 ());
2298         }
2299 
2300     private:
2301         matrix_closure_type data_;
2302         size_type lower_;
2303         size_type upper_;
2304         typedef const value_type const_value_type;
2305         static const_value_type zero_;
2306     };
2307 
2308     // Specialization for temporary_traits
2309     template <class M>
2310     struct vector_temporary_traits< banded_adaptor<M> >
2311     : vector_temporary_traits< M > {} ;
2312     template <class M>
2313     struct vector_temporary_traits< const banded_adaptor<M> >
2314     : vector_temporary_traits< M > {} ;
2315 
2316     template <class M>
2317     struct matrix_temporary_traits< banded_adaptor<M> >
2318     : matrix_temporary_traits< M > {} ;
2319     template <class M>
2320     struct matrix_temporary_traits< const banded_adaptor<M> >
2321     : matrix_temporary_traits< M > {} ;
2322 
2323 
2324     template<class M>
2325     typename banded_adaptor<M>::const_value_type banded_adaptor<M>::zero_ = value_type/*zero*/();
2326 
2327     /** \brief A diagonal matrix adaptator: convert a any matrix into a diagonal matrix expression
2328      *
2329      * For a \f$(m\times m)\f$-dimensional matrix, the \c diagonal_adaptor will provide a diagonal matrix
2330      * with \f$0 \leq i < m\f$ and \f$0 \leq j < m\f$, if \f$i\neq j\f$ then \f$b_{i,j}=0\f$.
2331      *
2332      * Storage and location are based on those of the underlying matrix. This is important because
2333      * a \c diagonal_adaptor does not copy the matrix data to a new place. Therefore, modifying values
2334      * in a \c diagonal_adaptor matrix will also modify the underlying matrix too.
2335      *
2336      * \tparam M the type of matrix used to generate the diagonal matrix
2337      */
2338 
2339     template<class M>
2340     class diagonal_adaptor:
2341         public banded_adaptor<M> {
2342     public:
2343         typedef M matrix_type;
2344         typedef banded_adaptor<M> adaptor_type;
2345 
2346         // Construction and destruction
2347         BOOST_UBLAS_INLINE
diagonal_adaptor()2348         diagonal_adaptor ():
2349             adaptor_type () {}
2350         BOOST_UBLAS_INLINE
diagonal_adaptor(matrix_type & data)2351         diagonal_adaptor (matrix_type &data):
2352             adaptor_type (data) {}
2353         BOOST_UBLAS_INLINE
~diagonal_adaptor()2354         ~diagonal_adaptor () {}
2355 
2356         // Assignment
2357         BOOST_UBLAS_INLINE
operator =(const diagonal_adaptor & m)2358         diagonal_adaptor &operator = (const diagonal_adaptor &m) {
2359             adaptor_type::operator = (m);
2360             return *this;
2361         }
2362         template<class AE>
2363         BOOST_UBLAS_INLINE
operator =(const matrix_expression<AE> & ae)2364         diagonal_adaptor &operator = (const matrix_expression<AE> &ae) {
2365             adaptor_type::operator = (ae);
2366             return *this;
2367         }
2368     };
2369 
2370 }}}
2371 
2372 #endif
2373