New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

matrix_proxy.hpp in vendors/XIOS/current/extern/boost/include/boost/numeric/ublas – NEMO

Context Navigation

source: vendors/XIOS/current/extern/boost/include/boost/numeric/ublas/matrix_proxy.hpp @ 3408

Last change on this file since 3408 was 3408, checked in by rblod, 12 years ago
importing initial XIOS vendor drop
Property svn:keywords set to `Id`
File size: 187.0 KB

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

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

Download in other formats: