Context Navigation

source: vendor/nemo/current/NEMOGCM/EXTERNAL/XIOS/extern/boost/include/boost/numeric/ublas/matrix_expression.hpp @ 44

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

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

Download in other formats: