Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

matrix_proxy.hpp @ 857

Revision 857, 186.1 KB checked in by igarcia, 18 years ago (diff)

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

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: NonGTP/Boost/boost/numeric/ublas/matrix_proxy.hpp @ 857

Download in other formats: