Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

matrix_assign.hpp @ 857

Revision 857, 78.1 KB checked in by igarcia, 19 years ago (diff)

Line
1	//
2	// Copyright (c) 2000-2002
3	// Joerg Walter, Mathias Koch
4	//
5	// Permission to use, copy, modify, distribute and sell this software
6	// and its documentation for any purpose is hereby granted without fee,
7	// provided that the above copyright notice appear in all copies and
8	// that both that copyright notice and this permission notice appear
9	// in supporting documentation. The authors make no representations
10	// about the suitability of this software for any purpose.
11	// It is provided "as is" without express or implied warranty.
12	//
13	// The authors gratefully acknowledge the support of
14	// GeNeSys mbH & Co. KG in producing this work.
15	//
16
17	#ifndef _BOOST_UBLAS_MATRIX_ASSIGN_
18	#define _BOOST_UBLAS_MATRIX_ASSIGN_
19
20	// Required for make_conformant storage
21	#include <vector>
22
23	// Iterators based on ideas of Jeremy Siek
24
25	namespace boost { namespace numeric { namespace ublas {
26	namespace detail {
27
28	// Weak equality check - useful to compare equality two arbitary matrix expression results.
29	// Since the actual expressions are unknown, we check for and arbitary error bound
30	// on the relative error.
31	// For a linear expression the infinity norm makes sense as we do not know how the elements will be
32	// combined in the expression. False positive results are inevitable for arbirary expressions!
33	template<class E1, class E2, class S>
34	BOOST_UBLAS_INLINE
35	bool equals (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2, S epsilon, S min_norm) {
36	return norm_inf (e1 - e2) < epsilon *
37	std::max<S> (std::max<S> (norm_inf (e1), norm_inf (e2)), min_norm);
38	}
39
40	template<class E1, class E2>
41	BOOST_UBLAS_INLINE
42	bool expression_type_check (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2) {
43	typedef typename type_traits<typename promote_traits<typename E1::value_type,
44	typename E2::value_type>::promote_type>::real_type real_type;
45	return equals (e1, e2, BOOST_UBLAS_TYPE_CHECK_EPSILON, BOOST_UBLAS_TYPE_CHECK_MIN);
46	}
47
48
49	template<class M, class E, class R>
50	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
51	void make_conformant (M &m, const matrix_expression<E> &e, row_major_tag, R) {
52	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
53	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
54	typedef R conformant_restrict_type;
55	typedef typename M::size_type size_type;
56	typedef typename M::difference_type difference_type;
57	typedef typename M::value_type value_type;
58	// FIXME unbounded_array with push_back maybe better
59	std::vector<std::pair<size_type, size_type> > index;
60	typename M::iterator1 it1 (m.begin1 ());
61	typename M::iterator1 it1_end (m.end1 ());
62	typename E::const_iterator1 it1e (e ().begin1 ());
63	typename E::const_iterator1 it1e_end (e ().end1 ());
64	while (it1 != it1_end && it1e != it1e_end) {
65	difference_type compare = it1.index1 () - it1e.index1 ();
66	if (compare == 0) {
67	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
68	typename M::iterator2 it2 (it1.begin ());
69	typename M::iterator2 it2_end (it1.end ());
70	typename E::const_iterator2 it2e (it1e.begin ());
71	typename E::const_iterator2 it2e_end (it1e.end ());
72	#else
73	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
74	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
75	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
76	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
77	#endif
78	if (it2 != it2_end && it2e != it2e_end) {
79	size_type it2_index = it2.index2 (), it2e_index = it2e.index2 ();
80	while (true) {
81	difference_type compare = it2_index - it2e_index;
82	if (compare == 0) {
83	++ it2, ++ it2e;
84	if (it2 != it2_end && it2e != it2e_end) {
85	it2_index = it2.index2 ();
86	it2e_index = it2e.index2 ();
87	} else
88	break;
89	} else if (compare < 0) {
90	increment (it2, it2_end, - compare);
91	if (it2 != it2_end)
92	it2_index = it2.index2 ();
93	else
94	break;
95	} else if (compare > 0) {
96	if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ()))
97	if (it2e != value_type/zero*/())
98	index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ()));
99	++ it2e;
100	if (it2e != it2e_end)
101	it2e_index = it2e.index2 ();
102	else
103	break;
104	}
105	}
106	}
107	while (it2e != it2e_end) {
108	if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ()))
109	if (it2e != value_type/zero*/())
110	index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ()));
111	++ it2e;
112	}
113	++ it1, ++ it1e;
114	} else if (compare < 0) {
115	increment (it1, it1_end, - compare);
116	} else if (compare > 0) {
117	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
118	typename E::const_iterator2 it2e (it1e.begin ());
119	typename E::const_iterator2 it2e_end (it1e.end ());
120	#else
121	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
122	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
123	#endif
124	while (it2e != it2e_end) {
125	if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ()))
126	if (it2e != value_type/zero*/())
127	index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ()));
128	++ it2e;
129	}
130	++ it1e;
131	}
132	}
133	while (it1e != it1e_end) {
134	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
135	typename E::const_iterator2 it2e (it1e.begin ());
136	typename E::const_iterator2 it2e_end (it1e.end ());
137	#else
138	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
139	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
140	#endif
141	while (it2e != it2e_end) {
142	if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ()))
143	if (it2e != value_type/zero*/())
144	index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ()));
145	++ it2e;
146	}
147	++ it1e;
148	}
149	// ISSUE proxies require insert_element
150	for (size_type k = 0; k < index.size (); ++ k)
151	m (index [k].first, index [k].second) = value_type/zero/();
152	}
153	template<class M, class E, class R>
154	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
155	void make_conformant (M &m, const matrix_expression<E> &e, column_major_tag, R) {
156	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
157	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
158	typedef R conformant_restrict_type;
159	typedef typename M::size_type size_type;
160	typedef typename M::difference_type difference_type;
161	typedef typename M::value_type value_type;
162	std::vector<std::pair<size_type, size_type> > index;
163	typename M::iterator2 it2 (m.begin2 ());
164	typename M::iterator2 it2_end (m.end2 ());
165	typename E::const_iterator2 it2e (e ().begin2 ());
166	typename E::const_iterator2 it2e_end (e ().end2 ());
167	while (it2 != it2_end && it2e != it2e_end) {
168	difference_type compare = it2.index2 () - it2e.index2 ();
169	if (compare == 0) {
170	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
171	typename M::iterator1 it1 (it2.begin ());
172	typename M::iterator1 it1_end (it2.end ());
173	typename E::const_iterator1 it1e (it2e.begin ());
174	typename E::const_iterator1 it1e_end (it2e.end ());
175	#else
176	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
177	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
178	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
179	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
180	#endif
181	if (it1 != it1_end && it1e != it1e_end) {
182	size_type it1_index = it1.index1 (), it1e_index = it1e.index1 ();
183	while (true) {
184	difference_type compare = it1_index - it1e_index;
185	if (compare == 0) {
186	++ it1, ++ it1e;
187	if (it1 != it1_end && it1e != it1e_end) {
188	it1_index = it1.index1 ();
189	it1e_index = it1e.index1 ();
190	} else
191	break;
192	} else if (compare < 0) {
193	increment (it1, it1_end, - compare);
194	if (it1 != it1_end)
195	it1_index = it1.index1 ();
196	else
197	break;
198	} else if (compare > 0) {
199	if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ()))
200	if (it1e != value_type/zero*/())
201	index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ()));
202	++ it1e;
203	if (it1e != it1e_end)
204	it1e_index = it1e.index1 ();
205	else
206	break;
207	}
208	}
209	}
210	while (it1e != it1e_end) {
211	if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ()))
212	if (it1e != value_type/zero*/())
213	index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ()));
214	++ it1e;
215	}
216	++ it2, ++ it2e;
217	} else if (compare < 0) {
218	increment (it2, it2_end, - compare);
219	} else if (compare > 0) {
220	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
221	typename E::const_iterator1 it1e (it2e.begin ());
222	typename E::const_iterator1 it1e_end (it2e.end ());
223	#else
224	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
225	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
226	#endif
227	while (it1e != it1e_end) {
228	if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ()))
229	if (it1e != value_type/zero*/())
230	index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ()));
231	++ it1e;
232	}
233	++ it2e;
234	}
235	}
236	while (it2e != it2e_end) {
237	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
238	typename E::const_iterator1 it1e (it2e.begin ());
239	typename E::const_iterator1 it1e_end (it2e.end ());
240	#else
241	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
242	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
243	#endif
244	while (it1e != it1e_end) {
245	if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ()))
246	if (it1e != value_type/zero*/())
247	index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ()));
248	++ it1e;
249	}
250	++ it2e;
251	}
252	// ISSUE proxies require insert_element
253	for (size_type k = 0; k < index.size (); ++ k)
254	m (index [k].first, index [k].second) = value_type/zero/();
255	}
256
257	}//namespace detail
258
259
260	// Explicitly iterating row major
261	template<template <class T1, class T2> class F, class M, class T>
262	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
263	void iterating_matrix_assign_scalar (M &m, const T &t, row_major_tag) {
264	typedef F<typename M::iterator2::reference, T> functor_type;
265	typedef typename M::difference_type difference_type;
266	difference_type size1 (m.size1 ());
267	difference_type size2 (m.size2 ());
268	typename M::iterator1 it1 (m.begin1 ());
269	BOOST_UBLAS_CHECK (size2 == 0 \|\| m.end1 () - it1 == size1, bad_size ());
270	while (-- size1 >= 0) {
271	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
272	typename M::iterator2 it2 (it1.begin ());
273	#else
274	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
275	#endif
276	BOOST_UBLAS_CHECK (it1.end () - it2 == size2, bad_size ());
277	difference_type temp_size2 (size2);
278	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
279	while (-- temp_size2 >= 0)
280	functor_type::apply (*it2, t), ++ it2;
281	#else
282	DD (temp_size2, 4, r, (functor_type::apply (*it2, t), ++ it2));
283	#endif
284	++ it1;
285	}
286	}
287	// Explicitly iterating column major
288	template<template <class T1, class T2> class F, class M, class T>
289	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
290	void iterating_matrix_assign_scalar (M &m, const T &t, column_major_tag) {
291	typedef F<typename M::iterator1::reference, T> functor_type;
292	typedef typename M::difference_type difference_type;
293	difference_type size2 (m.size2 ());
294	difference_type size1 (m.size1 ());
295	typename M::iterator2 it2 (m.begin2 ());
296	BOOST_UBLAS_CHECK (size1 == 0 \|\| m.end2 () - it2 == size2, bad_size ());
297	while (-- size2 >= 0) {
298	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
299	typename M::iterator1 it1 (it2.begin ());
300	#else
301	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
302	#endif
303	BOOST_UBLAS_CHECK (it2.end () - it1 == size1, bad_size ());
304	difference_type temp_size1 (size1);
305	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
306	while (-- temp_size1 >= 0)
307	functor_type::apply (*it1, t), ++ it1;
308	#else
309	DD (temp_size1, 4, r, (functor_type::apply (*it1, t), ++ it1));
310	#endif
311	++ it2;
312	}
313	}
314	// Explicitly indexing row major
315	template<template <class T1, class T2> class F, class M, class T>
316	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
317	void indexing_matrix_assign_scalar (M &m, const T &t, row_major_tag) {
318	typedef F<typename M::reference, T> functor_type;
319	typedef typename M::size_type size_type;
320	size_type size1 (m.size1 ());
321	size_type size2 (m.size2 ());
322	for (size_type i = 0; i < size1; ++ i) {
323	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
324	for (size_type j = 0; j < size2; ++ j)
325	functor_type::apply (m (i, j), t);
326	#else
327	size_type j (0);
328	DD (size2, 4, r, (functor_type::apply (m (i, j), t), ++ j));
329	#endif
330	}
331	}
332	// Explicitly indexing column major
333	template<template <class T1, class T2> class F, class M, class T>
334	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
335	void indexing_matrix_assign_scalar (M &m, const T &t, column_major_tag) {
336	typedef F<typename M::reference, T> functor_type;
337	typedef typename M::size_type size_type;
338	size_type size2 (m.size2 ());
339	size_type size1 (m.size1 ());
340	for (size_type j = 0; j < size2; ++ j) {
341	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
342	for (size_type i = 0; i < size1; ++ i)
343	functor_type::apply (m (i, j), t);
344	#else
345	size_type i (0);
346	DD (size1, 4, r, (functor_type::apply (m (i, j), t), ++ i));
347	#endif
348	}
349	}
350
351	// Dense (proxy) case
352	template<template <class T1, class T2> class F, class M, class T, class C>
353	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
354	void matrix_assign_scalar (M &m, const T &t, dense_proxy_tag, C) {
355	typedef C orientation_category;
356	#ifdef BOOST_UBLAS_USE_INDEXING
357	indexing_matrix_assign_scalar<F> (m, t, orientation_category ());
358	#elif BOOST_UBLAS_USE_ITERATING
359	iterating_matrix_assign_scalar<F> (m, t, orientation_category ());
360	#else
361	typedef typename M::size_type size_type;
362	size_type size1 (m.size1 ());
363	size_type size2 (m.size2 ());
364	if (size1 >= BOOST_UBLAS_ITERATOR_THRESHOLD &&
365	size2 >= BOOST_UBLAS_ITERATOR_THRESHOLD)
366	iterating_matrix_assign_scalar<F> (m, t, orientation_category ());
367	else
368	indexing_matrix_assign_scalar<F> (m, t, orientation_category ());
369	#endif
370	}
371	// Packed (proxy) row major case
372	template<template <class T1, class T2> class F, class M, class T>
373	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
374	void matrix_assign_scalar (M &m, const T &t, packed_proxy_tag, row_major_tag) {
375	typedef F<typename M::iterator2::reference, T> functor_type;
376	typedef typename M::difference_type difference_type;
377	typename M::iterator1 it1 (m.begin1 ());
378	difference_type size1 (m.end1 () - it1);
379	while (-- size1 >= 0) {
380	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
381	typename M::iterator2 it2 (it1.begin ());
382	difference_type size2 (it1.end () - it2);
383	#else
384	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
385	difference_type size2 (end (it1, iterator1_tag ()) - it2);
386	#endif
387	while (-- size2 >= 0)
388	functor_type::apply (*it2, t), ++ it2;
389	++ it1;
390	}
391	}
392	// Packed (proxy) column major case
393	template<template <class T1, class T2> class F, class M, class T>
394	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
395	void matrix_assign_scalar (M &m, const T &t, packed_proxy_tag, column_major_tag) {
396	typedef F<typename M::iterator1::reference, T> functor_type;
397	typedef typename M::difference_type difference_type;
398	typename M::iterator2 it2 (m.begin2 ());
399	difference_type size2 (m.end2 () - it2);
400	while (-- size2 >= 0) {
401	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
402	typename M::iterator1 it1 (it2.begin ());
403	difference_type size1 (it2.end () - it1);
404	#else
405	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
406	difference_type size1 (end (it2, iterator2_tag ()) - it1);
407	#endif
408	while (-- size1 >= 0)
409	functor_type::apply (*it1, t), ++ it1;
410	++ it2;
411	}
412	}
413	// Sparse (proxy) row major case
414	template<template <class T1, class T2> class F, class M, class T>
415	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
416	void matrix_assign_scalar (M &m, const T &t, sparse_proxy_tag, row_major_tag) {
417	typedef F<typename M::iterator2::reference, T> functor_type;
418	typename M::iterator1 it1 (m.begin1 ());
419	typename M::iterator1 it1_end (m.end1 ());
420	while (it1 != it1_end) {
421	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
422	typename M::iterator2 it2 (it1.begin ());
423	typename M::iterator2 it2_end (it1.end ());
424	#else
425	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
426	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
427	#endif
428	while (it2 != it2_end)
429	functor_type::apply (*it2, t), ++ it2;
430	++ it1;
431	}
432	}
433	// Sparse (proxy) column major case
434	template<template <class T1, class T2> class F, class M, class T>
435	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
436	void matrix_assign_scalar (M &m, const T &t, sparse_proxy_tag, column_major_tag) {
437	typedef F<typename M::iterator1::reference, T> functor_type;
438	typename M::iterator2 it2 (m.begin2 ());
439	typename M::iterator2 it2_end (m.end2 ());
440	while (it2 != it2_end) {
441	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
442	typename M::iterator1 it1 (it2.begin ());
443	typename M::iterator1 it1_end (it2.end ());
444	#else
445	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
446	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
447	#endif
448	while (it1 != it1_end)
449	functor_type::apply (*it1, t), ++ it1;
450	++ it2;
451	}
452	}
453
454	// Dispatcher
455	template<template <class T1, class T2> class F, class M, class T>
456	BOOST_UBLAS_INLINE
457	void matrix_assign_scalar (M &m, const T &t) {
458	typedef typename M::storage_category storage_category;
459	typedef typename M::orientation_category orientation_category;
460	matrix_assign_scalar<F> (m, t, storage_category (), orientation_category ());
461	}
462
463	template<class SC, bool COMPUTED, class RI1, class RI2>
464	struct matrix_assign_traits {
465	typedef SC storage_category;
466	};
467
468	template<bool COMPUTED>
469	struct matrix_assign_traits<dense_tag, COMPUTED, packed_random_access_iterator_tag, packed_random_access_iterator_tag> {
470	typedef packed_tag storage_category;
471	};
472	template<>
473	struct matrix_assign_traits<dense_tag, false, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
474	typedef sparse_tag storage_category;
475	};
476	template<>
477	struct matrix_assign_traits<dense_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
478	typedef sparse_proxy_tag storage_category;
479	};
480
481	template<bool COMPUTED>
482	struct matrix_assign_traits<dense_proxy_tag, COMPUTED, packed_random_access_iterator_tag, packed_random_access_iterator_tag> {
483	typedef packed_proxy_tag storage_category;
484	};
485	template<bool COMPUTED>
486	struct matrix_assign_traits<dense_proxy_tag, COMPUTED, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
487	typedef sparse_proxy_tag storage_category;
488	};
489
490	template<>
491	struct matrix_assign_traits<packed_tag, false, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
492	typedef sparse_tag storage_category;
493	};
494	template<>
495	struct matrix_assign_traits<packed_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
496	typedef sparse_proxy_tag storage_category;
497	};
498
499	template<bool COMPUTED>
500	struct matrix_assign_traits<packed_proxy_tag, COMPUTED, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
501	typedef sparse_proxy_tag storage_category;
502	};
503
504	template<>
505	struct matrix_assign_traits<sparse_tag, true, dense_random_access_iterator_tag, dense_random_access_iterator_tag> {
506	typedef sparse_proxy_tag storage_category;
507	};
508	template<>
509	struct matrix_assign_traits<sparse_tag, true, packed_random_access_iterator_tag, packed_random_access_iterator_tag> {
510	typedef sparse_proxy_tag storage_category;
511	};
512	template<>
513	struct matrix_assign_traits<sparse_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
514	typedef sparse_proxy_tag storage_category;
515	};
516
517	// Explicitly iterating row major
518	template<template <class T1, class T2> class F, class M, class E>
519	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
520	void iterating_matrix_assign (M &m, const matrix_expression<E> &e, row_major_tag) {
521	typedef F<typename M::iterator2::reference, typename E::value_type> functor_type;
522	typedef typename M::difference_type difference_type;
523	difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ()));
524	difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ()));
525	typename M::iterator1 it1 (m.begin1 ());
526	BOOST_UBLAS_CHECK (size2 == 0 \|\| m.end1 () - it1 == size1, bad_size ());
527	typename E::const_iterator1 it1e (e ().begin1 ());
528	BOOST_UBLAS_CHECK (size2 == 0 \|\| e ().end1 () - it1e == size1, bad_size ());
529	while (-- size1 >= 0) {
530	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
531	typename M::iterator2 it2 (it1.begin ());
532	typename E::const_iterator2 it2e (it1e.begin ());
533	#else
534	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
535	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
536	#endif
537	BOOST_UBLAS_CHECK (it1.end () - it2 == size2, bad_size ());
538	BOOST_UBLAS_CHECK (it1e.end () - it2e == size2, bad_size ());
539	difference_type temp_size2 (size2);
540	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
541	while (-- temp_size2 >= 0)
542	functor_type::apply (it2, it2e), ++ it2, ++ it2e;
543	#else
544	DD (temp_size2, 2, r, (functor_type::apply (it2, it2e), ++ it2, ++ it2e));
545	#endif
546	++ it1, ++ it1e;
547	}
548	}
549	// Explicitly iterating column major
550	template<template <class T1, class T2> class F, class M, class E>
551	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
552	void iterating_matrix_assign (M &m, const matrix_expression<E> &e, column_major_tag) {
553	typedef F<typename M::iterator1::reference, typename E::value_type> functor_type;
554	typedef typename M::difference_type difference_type;
555	difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ()));
556	difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ()));
557	typename M::iterator2 it2 (m.begin2 ());
558	BOOST_UBLAS_CHECK (size1 == 0 \|\| m.end2 () - it2 == size2, bad_size ());
559	typename E::const_iterator2 it2e (e ().begin2 ());
560	BOOST_UBLAS_CHECK (size1 == 0 \|\| e ().end2 () - it2e == size2, bad_size ());
561	while (-- size2 >= 0) {
562	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
563	typename M::iterator1 it1 (it2.begin ());
564	typename E::const_iterator1 it1e (it2e.begin ());
565	#else
566	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
567	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
568	#endif
569	BOOST_UBLAS_CHECK (it2.end () - it1 == size1, bad_size ());
570	BOOST_UBLAS_CHECK (it2e.end () - it1e == size1, bad_size ());
571	difference_type temp_size1 (size1);
572	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
573	while (-- temp_size1 >= 0)
574	functor_type::apply (it1, it1e), ++ it1, ++ it1e;
575	#else
576	DD (temp_size1, 2, r, (functor_type::apply (it1, it1e), ++ it1, ++ it1e));
577	#endif
578	++ it2, ++ it2e;
579	}
580	}
581	// Explicitly indexing row major
582	template<template <class T1, class T2> class F, class M, class E>
583	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
584	void indexing_matrix_assign (M &m, const matrix_expression<E> &e, row_major_tag) {
585	typedef F<typename M::reference, typename E::value_type> functor_type;
586	typedef typename M::size_type size_type;
587	size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ()));
588	size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ()));
589	for (size_type i = 0; i < size1; ++ i) {
590	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
591	for (size_type j = 0; j < size2; ++ j)
592	functor_type::apply (m (i, j), e () (i, j));
593	#else
594	size_type j (0);
595	DD (size2, 2, r, (functor_type::apply (m (i, j), e () (i, j)), ++ j));
596	#endif
597	}
598	}
599	// Explicitly indexing column major
600	template<template <class T1, class T2> class F, class M, class E>
601	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
602	void indexing_matrix_assign (M &m, const matrix_expression<E> &e, column_major_tag) {
603	typedef F<typename M::reference, typename E::value_type> functor_type;
604	typedef typename M::size_type size_type;
605	size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ()));
606	size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ()));
607	for (size_type j = 0; j < size2; ++ j) {
608	#ifndef BOOST_UBLAS_USE_DUFF_DEVICE
609	for (size_type i = 0; i < size1; ++ i)
610	functor_type::apply (m (i, j), e () (i, j));
611	#else
612	size_type i (0);
613	DD (size1, 2, r, (functor_type::apply (m (i, j), e () (i, j)), ++ i));
614	#endif
615	}
616	}
617
618	// Dense (proxy) case
619	template<template <class T1, class T2> class F, class R, class M, class E, class C>
620	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
621	void matrix_assign (M &m, const matrix_expression<E> &e, dense_proxy_tag, C) {
622	// R unnecessary, make_conformant not required
623	typedef C orientation_category;
624	#ifdef BOOST_UBLAS_USE_INDEXING
625	indexing_matrix_assign<F> (m, e, orientation_category ());
626	#elif BOOST_UBLAS_USE_ITERATING
627	iterating_matrix_assign<F> (m, e, orientation_category ());
628	#else
629	typedef typename M::difference_type difference_type;
630	size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ()));
631	size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ()));
632	if (size1 >= BOOST_UBLAS_ITERATOR_THRESHOLD &&
633	size2 >= BOOST_UBLAS_ITERATOR_THRESHOLD)
634	iterating_matrix_assign<F> (m, e, orientation_category ());
635	else
636	indexing_matrix_assign<F> (m, e, orientation_category ());
637	#endif
638	}
639	// Packed (proxy) row major case
640	template<template <class T1, class T2> class F, class R, class M, class E>
641	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
642	void matrix_assign (M &m, const matrix_expression<E> &e, packed_proxy_tag, row_major_tag) {
643	typedef F<typename M::iterator2::reference, typename E::value_type> functor_type;
644	// R unnecessary, make_conformant not required
645	typedef typename M::difference_type difference_type;
646	typedef typename M::value_type value_type;
647	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
648	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
649	#if BOOST_UBLAS_TYPE_CHECK
650	matrix<value_type, row_major> cm (m.size1 (), m.size2 ());
651	indexing_matrix_assign<scalar_assign> (cm, m, row_major_tag ());
652	indexing_matrix_assign<F> (cm, e, row_major_tag ());
653	#endif
654	typename M::iterator1 it1 (m.begin1 ());
655	typename M::iterator1 it1_end (m.end1 ());
656	typename E::const_iterator1 it1e (e ().begin1 ());
657	typename E::const_iterator1 it1e_end (e ().end1 ());
658	difference_type it1_size (it1_end - it1);
659	difference_type it1e_size (it1e_end - it1e);
660	difference_type diff1 (0);
661	if (it1_size > 0 && it1e_size > 0)
662	diff1 = it1.index1 () - it1e.index1 ();
663	if (diff1 != 0) {
664	difference_type size1 = (std::min) (diff1, it1e_size);
665	if (size1 > 0) {
666	it1e += size1;
667	it1e_size -= size1;
668	diff1 -= size1;
669	}
670	size1 = (std::min) (- diff1, it1_size);
671	if (size1 > 0) {
672	it1_size -= size1;
673	if (!functor_type::computed) {
674	while (-- size1 >= 0) { // zeroing
675	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
676	typename M::iterator2 it2 (it1.begin ());
677	typename M::iterator2 it2_end (it1.end ());
678	#else
679	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
680	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
681	#endif
682	difference_type size2 (it2_end - it2);
683	while (-- size2 >= 0)
684	functor_type::apply (it2, value_type/zero*/()), ++ it2;
685	++ it1;
686	}
687	} else {
688	it1 += size1;
689	}
690	diff1 += size1;
691	}
692	}
693	difference_type size1 ((std::min) (it1_size, it1e_size));
694	it1_size -= size1;
695	it1e_size -= size1;
696	while (-- size1 >= 0) {
697	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
698	typename M::iterator2 it2 (it1.begin ());
699	typename M::iterator2 it2_end (it1.end ());
700	typename E::const_iterator2 it2e (it1e.begin ());
701	typename E::const_iterator2 it2e_end (it1e.end ());
702	#else
703	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
704	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
705	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
706	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
707	#endif
708	difference_type it2_size (it2_end - it2);
709	difference_type it2e_size (it2e_end - it2e);
710	difference_type diff2 (0);
711	if (it2_size > 0 && it2e_size > 0) {
712	diff2 = it2.index2 () - it2e.index2 ();
713	difference_type size2 = (std::min) (diff2, it2e_size);
714	if (size2 > 0) {
715	it2e += size2;
716	it2e_size -= size2;
717	diff2 -= size2;
718	}
719	size2 = (std::min) (- diff2, it2_size);
720	if (size2 > 0) {
721	it2_size -= size2;
722	if (!functor_type::computed) {
723	while (-- size2 >= 0) // zeroing
724	functor_type::apply (it2, value_type/zero*/()), ++ it2;
725	} else {
726	it2 += size2;
727	}
728	diff2 += size2;
729	}
730	}
731	difference_type size2 ((std::min) (it2_size, it2e_size));
732	it2_size -= size2;
733	it2e_size -= size2;
734	while (-- size2 >= 0)
735	functor_type::apply (it2, it2e), ++ it2, ++ it2e;
736	size2 = it2_size;
737	if (!functor_type::computed) {
738	while (-- size2 >= 0) // zeroing
739	functor_type::apply (it2, value_type/zero*/()), ++ it2;
740	} else {
741	it2 += size2;
742	}
743	++ it1, ++ it1e;
744	}
745	size1 = it1_size;
746	if (!functor_type::computed) {
747	while (-- size1 >= 0) { // zeroing
748	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
749	typename M::iterator2 it2 (it1.begin ());
750	typename M::iterator2 it2_end (it1.end ());
751	#else
752	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
753	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
754	#endif
755	difference_type size2 (it2_end - it2);
756	while (-- size2 >= 0)
757	functor_type::apply (it2, value_type/zero*/()), ++ it2;
758	++ it1;
759	}
760	} else {
761	it1 += size1;
762	}
763	#if BOOST_UBLAS_TYPE_CHECK
764	if (! disable_type_check<bool>::value)
765	BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ());
766	#endif
767	}
768	// Packed (proxy) column major case
769	template<template <class T1, class T2> class F, class R, class M, class E>
770	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
771	void matrix_assign (M &m, const matrix_expression<E> &e, packed_proxy_tag, column_major_tag) {
772	typedef F<typename M::iterator1::reference, typename E::value_type> functor_type;
773	// R unnecessary, make_conformant not required
774	typedef typename M::difference_type difference_type;
775	typedef typename M::value_type value_type;
776	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
777	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
778	#if BOOST_UBLAS_TYPE_CHECK
779	matrix<value_type, column_major> cm (m.size1 (), m.size2 ());
780	indexing_matrix_assign<scalar_assign> (cm, m, column_major_tag ());
781	indexing_matrix_assign<F> (cm, e, column_major_tag ());
782	#endif
783	typename M::iterator2 it2 (m.begin2 ());
784	typename M::iterator2 it2_end (m.end2 ());
785	typename E::const_iterator2 it2e (e ().begin2 ());
786	typename E::const_iterator2 it2e_end (e ().end2 ());
787	difference_type it2_size (it2_end - it2);
788	difference_type it2e_size (it2e_end - it2e);
789	difference_type diff2 (0);
790	if (it2_size > 0 && it2e_size > 0)
791	diff2 = it2.index2 () - it2e.index2 ();
792	if (diff2 != 0) {
793	difference_type size2 = (std::min) (diff2, it2e_size);
794	if (size2 > 0) {
795	it2e += size2;
796	it2e_size -= size2;
797	diff2 -= size2;
798	}
799	size2 = (std::min) (- diff2, it2_size);
800	if (size2 > 0) {
801	it2_size -= size2;
802	if (!functor_type::computed) {
803	while (-- size2 >= 0) { // zeroing
804	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
805	typename M::iterator1 it1 (it2.begin ());
806	typename M::iterator1 it1_end (it2.end ());
807	#else
808	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
809	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
810	#endif
811	difference_type size1 (it1_end - it1);
812	while (-- size1 >= 0)
813	functor_type::apply (it1, value_type/zero*/()), ++ it1;
814	++ it2;
815	}
816	} else {
817	it2 += size2;
818	}
819	diff2 += size2;
820	}
821	}
822	difference_type size2 ((std::min) (it2_size, it2e_size));
823	it2_size -= size2;
824	it2e_size -= size2;
825	while (-- size2 >= 0) {
826	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
827	typename M::iterator1 it1 (it2.begin ());
828	typename M::iterator1 it1_end (it2.end ());
829	typename E::const_iterator1 it1e (it2e.begin ());
830	typename E::const_iterator1 it1e_end (it2e.end ());
831	#else
832	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
833	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
834	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
835	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
836	#endif
837	difference_type it1_size (it1_end - it1);
838	difference_type it1e_size (it1e_end - it1e);
839	difference_type diff1 (0);
840	if (it1_size > 0 && it1e_size > 0) {
841	diff1 = it1.index1 () - it1e.index1 ();
842	difference_type size1 = (std::min) (diff1, it1e_size);
843	if (size1 > 0) {
844	it1e += size1;
845	it1e_size -= size1;
846	diff1 -= size1;
847	}
848	size1 = (std::min) (- diff1, it1_size);
849	if (size1 > 0) {
850	it1_size -= size1;
851	if (!functor_type::computed) {
852	while (-- size1 >= 0) // zeroing
853	functor_type::apply (it1, value_type/zero*/()), ++ it1;
854	} else {
855	it1 += size1;
856	}
857	diff1 += size1;
858	}
859	}
860	difference_type size1 ((std::min) (it1_size, it1e_size));
861	it1_size -= size1;
862	it1e_size -= size1;
863	while (-- size1 >= 0)
864	functor_type::apply (it1, it1e), ++ it1, ++ it1e;
865	size1 = it1_size;
866	if (!functor_type::computed) {
867	while (-- size1 >= 0) // zeroing
868	functor_type::apply (it1, value_type/zero*/()), ++ it1;
869	} else {
870	it1 += size1;
871	}
872	++ it2, ++ it2e;
873	}
874	size2 = it2_size;
875	if (!functor_type::computed) {
876	while (-- size2 >= 0) { // zeroing
877	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
878	typename M::iterator1 it1 (it2.begin ());
879	typename M::iterator1 it1_end (it2.end ());
880	#else
881	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
882	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
883	#endif
884	difference_type size1 (it1_end - it1);
885	while (-- size1 >= 0)
886	functor_type::apply (it1, value_type/zero*/()), ++ it1;
887	++ it2;
888	}
889	} else {
890	it2 += size2;
891	}
892	#if BOOST_UBLAS_TYPE_CHECK
893	if (! disable_type_check<bool>::value)
894	BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ());
895	#endif
896	}
897	// Sparse row major case
898	template<template <class T1, class T2> class F, class R, class M, class E>
899	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
900	void matrix_assign (M &m, const matrix_expression<E> &e, sparse_tag, row_major_tag) {
901	typedef F<typename M::iterator2::reference, typename E::value_type> functor_type;
902	// R unnecessary, make_conformant not required
903	BOOST_STATIC_ASSERT ((!functor_type::computed));
904	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
905	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
906	typedef typename M::value_type value_type;
907	// Sparse type has no numeric constraints to check
908
909	m.clear ();
910	typename E::const_iterator1 it1e (e ().begin1 ());
911	typename E::const_iterator1 it1e_end (e ().end1 ());
912	while (it1e != it1e_end) {
913	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
914	typename E::const_iterator2 it2e (it1e.begin ());
915	typename E::const_iterator2 it2e_end (it1e.end ());
916	#else
917	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
918	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
919	#endif
920	while (it2e != it2e_end) {
921	value_type t (*it2e);
922	if (t != value_type/zero/())
923	m.insert_element (it2e.index1 (), it2e.index2 (), t);
924	++ it2e;
925	}
926	++ it1e;
927	}
928	}
929	// Sparse column major case
930	template<template <class T1, class T2> class F, class R, class M, class E>
931	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
932	void matrix_assign (M &m, const matrix_expression<E> &e, sparse_tag, column_major_tag) {
933	typedef F<typename M::iterator1::reference, typename E::value_type> functor_type;
934	// R unnecessary, make_conformant not required
935	BOOST_STATIC_ASSERT ((!functor_type::computed));
936	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
937	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
938	typedef typename M::value_type value_type;
939	// Sparse type has no numeric constraints to check
940
941	m.clear ();
942	typename E::const_iterator2 it2e (e ().begin2 ());
943	typename E::const_iterator2 it2e_end (e ().end2 ());
944	while (it2e != it2e_end) {
945	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
946	typename E::const_iterator1 it1e (it2e.begin ());
947	typename E::const_iterator1 it1e_end (it2e.end ());
948	#else
949	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
950	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
951	#endif
952	while (it1e != it1e_end) {
953	value_type t (*it1e);
954	if (t != value_type/zero/())
955	m.insert_element (it1e.index1 (), it1e.index2 (), t);
956	++ it1e;
957	}
958	++ it2e;
959	}
960	}
961	// Sparse proxy or functional row major case
962	template<template <class T1, class T2> class F, class R, class M, class E>
963	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
964	void matrix_assign (M &m, const matrix_expression<E> &e, sparse_proxy_tag, row_major_tag) {
965	typedef F<typename M::iterator2::reference, typename E::value_type> functor_type;
966	typedef R conformant_restrict_type;
967	typedef typename M::size_type size_type;
968	typedef typename M::difference_type difference_type;
969	typedef typename M::value_type value_type;
970	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
971	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
972	#if BOOST_UBLAS_TYPE_CHECK
973	matrix<value_type, row_major> cm (m.size1 (), m.size2 ());
974	indexing_matrix_assign<scalar_assign> (cm, m, row_major_tag ());
975	indexing_matrix_assign<F> (cm, e, row_major_tag ());
976	#endif
977	detail::make_conformant (m, e, row_major_tag (), conformant_restrict_type ());
978
979	typename M::iterator1 it1 (m.begin1 ());
980	typename M::iterator1 it1_end (m.end1 ());
981	typename E::const_iterator1 it1e (e ().begin1 ());
982	typename E::const_iterator1 it1e_end (e ().end1 ());
983	while (it1 != it1_end && it1e != it1e_end) {
984	difference_type compare = it1.index1 () - it1e.index1 ();
985	if (compare == 0) {
986	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
987	typename M::iterator2 it2 (it1.begin ());
988	typename M::iterator2 it2_end (it1.end ());
989	typename E::const_iterator2 it2e (it1e.begin ());
990	typename E::const_iterator2 it2e_end (it1e.end ());
991	#else
992	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
993	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
994	typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ()));
995	typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ()));
996	#endif
997	if (it2 != it2_end && it2e != it2e_end) {
998	size_type it2_index = it2.index2 (), it2e_index = it2e.index2 ();
999	while (true) {
1000	difference_type compare = it2_index - it2e_index;
1001	if (compare == 0) {
1002	functor_type::apply (it2, it2e);
1003	++ it2, ++ it2e;
1004	if (it2 != it2_end && it2e != it2e_end) {
1005	it2_index = it2.index2 ();
1006	it2e_index = it2e.index2 ();
1007	} else
1008	break;
1009	} else if (compare < 0) {
1010	if (!functor_type::computed) {
1011	functor_type::apply (it2, value_type/zero*/());
1012	++ it2;
1013	} else
1014	increment (it2, it2_end, - compare);
1015	if (it2 != it2_end)
1016	it2_index = it2.index2 ();
1017	else
1018	break;
1019	} else if (compare > 0) {
1020	increment (it2e, it2e_end, compare);
1021	if (it2e != it2e_end)
1022	it2e_index = it2e.index2 ();
1023	else
1024	break;
1025	}
1026	}
1027	}
1028	if (!functor_type::computed) {
1029	while (it2 != it2_end) { // zeroing
1030	functor_type::apply (it2, value_type/zero*/());
1031	++ it2;
1032	}
1033	} else {
1034	it2 = it2_end;
1035	}
1036	++ it1, ++ it1e;
1037	} else if (compare < 0) {
1038	if (!functor_type::computed) {
1039	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1040	typename M::iterator2 it2 (it1.begin ());
1041	typename M::iterator2 it2_end (it1.end ());
1042	#else
1043	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1044	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
1045	#endif
1046	while (it2 != it2_end) { // zeroing
1047	functor_type::apply (it2, value_type/zero*/());
1048	++ it2;
1049	}
1050	++ it1;
1051	} else {
1052	increment (it1, it1_end, - compare);
1053	}
1054	} else if (compare > 0) {
1055	increment (it1e, it1e_end, compare);
1056	}
1057	}
1058	if (!functor_type::computed) {
1059	while (it1 != it1_end) {
1060	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1061	typename M::iterator2 it2 (it1.begin ());
1062	typename M::iterator2 it2_end (it1.end ());
1063	#else
1064	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1065	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
1066	#endif
1067	while (it2 != it2_end) { // zeroing
1068	functor_type::apply (it2, value_type/zero*/());
1069	++ it2;
1070	}
1071	++ it1;
1072	}
1073	} else {
1074	it1 = it1_end;
1075	}
1076	#if BOOST_UBLAS_TYPE_CHECK
1077	if (! disable_type_check<bool>::value)
1078	BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ());
1079	#endif
1080	}
1081	// Sparse proxy or functional column major case
1082	template<template <class T1, class T2> class F, class R, class M, class E>
1083	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1084	void matrix_assign (M &m, const matrix_expression<E> &e, sparse_proxy_tag, column_major_tag) {
1085	typedef F<typename M::iterator1::reference, typename E::value_type> functor_type;
1086	typedef R conformant_restrict_type;
1087	typedef typename M::size_type size_type;
1088	typedef typename M::difference_type difference_type;
1089	typedef typename M::value_type value_type;
1090	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
1091	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
1092	#if BOOST_UBLAS_TYPE_CHECK
1093	matrix<value_type, column_major> cm (m.size1 (), m.size2 ());
1094	indexing_matrix_assign<scalar_assign> (cm, m, column_major_tag ());
1095	indexing_matrix_assign<F> (cm, e, column_major_tag ());
1096	#endif
1097	detail::make_conformant (m, e, column_major_tag (), conformant_restrict_type ());
1098
1099	typename M::iterator2 it2 (m.begin2 ());
1100	typename M::iterator2 it2_end (m.end2 ());
1101	typename E::const_iterator2 it2e (e ().begin2 ());
1102	typename E::const_iterator2 it2e_end (e ().end2 ());
1103	while (it2 != it2_end && it2e != it2e_end) {
1104	difference_type compare = it2.index2 () - it2e.index2 ();
1105	if (compare == 0) {
1106	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1107	typename M::iterator1 it1 (it2.begin ());
1108	typename M::iterator1 it1_end (it2.end ());
1109	typename E::const_iterator1 it1e (it2e.begin ());
1110	typename E::const_iterator1 it1e_end (it2e.end ());
1111	#else
1112	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1113	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1114	typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ()));
1115	typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ()));
1116	#endif
1117	if (it1 != it1_end && it1e != it1e_end) {
1118	size_type it1_index = it1.index1 (), it1e_index = it1e.index1 ();
1119	while (true) {
1120	difference_type compare = it1_index - it1e_index;
1121	if (compare == 0) {
1122	functor_type::apply (it1, it1e);
1123	++ it1, ++ it1e;
1124	if (it1 != it1_end && it1e != it1e_end) {
1125	it1_index = it1.index1 ();
1126	it1e_index = it1e.index1 ();
1127	} else
1128	break;
1129	} else if (compare < 0) {
1130	if (!functor_type::computed) {
1131	functor_type::apply (it1, value_type/zero*/()); // zeroing
1132	++ it1;
1133	} else
1134	increment (it1, it1_end, - compare);
1135	if (it1 != it1_end)
1136	it1_index = it1.index1 ();
1137	else
1138	break;
1139	} else if (compare > 0) {
1140	increment (it1e, it1e_end, compare);
1141	if (it1e != it1e_end)
1142	it1e_index = it1e.index1 ();
1143	else
1144	break;
1145	}
1146	}
1147	}
1148	if (!functor_type::computed) {
1149	while (it1 != it1_end) { // zeroing
1150	functor_type::apply (it1, value_type/zero*/());
1151	++ it1;
1152	}
1153	} else {
1154	it1 = it1_end;
1155	}
1156	++ it2, ++ it2e;
1157	} else if (compare < 0) {
1158	if (!functor_type::computed) {
1159	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1160	typename M::iterator1 it1 (it2.begin ());
1161	typename M::iterator1 it1_end (it2.end ());
1162	#else
1163	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1164	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1165	#endif
1166	while (it1 != it1_end) { // zeroing
1167	functor_type::apply (it1, value_type/zero*/());
1168	++ it1;
1169	}
1170	++ it2;
1171	} else {
1172	increment (it2, it2_end, - compare);
1173	}
1174	} else if (compare > 0) {
1175	increment (it2e, it2e_end, compare);
1176	}
1177	}
1178	if (!functor_type::computed) {
1179	while (it2 != it2_end) {
1180	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1181	typename M::iterator1 it1 (it2.begin ());
1182	typename M::iterator1 it1_end (it2.end ());
1183	#else
1184	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1185	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1186	#endif
1187	while (it1 != it1_end) { // zeroing
1188	functor_type::apply (it1, value_type/zero*/());
1189	++ it1;
1190	}
1191	++ it2;
1192	}
1193	} else {
1194	it2 = it2_end;
1195	}
1196	#if BOOST_UBLAS_TYPE_CHECK
1197	if (! disable_type_check<bool>::value)
1198	BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ());
1199	#endif
1200	}
1201
1202	// Dispatcher
1203	template<template <class T1, class T2> class F, class M, class E>
1204	BOOST_UBLAS_INLINE
1205	void matrix_assign (M &m, const matrix_expression<E> &e) {
1206	typedef typename matrix_assign_traits<typename M::storage_category,
1207	F<typename M::reference, typename E::value_type>::computed,
1208	typename E::const_iterator1::iterator_category,
1209	typename E::const_iterator2::iterator_category>::storage_category storage_category;
1210	// give preference to matrix M's orientation if known
1211	typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>,
1212	typename E::orientation_category ,
1213	typename M::orientation_category >::type orientation_category;
1214	typedef basic_full<typename M::size_type> unrestricted;
1215	matrix_assign<F, unrestricted> (m, e, storage_category (), orientation_category ());
1216	}
1217	template<template <class T1, class T2> class F, class R, class M, class E>
1218	BOOST_UBLAS_INLINE
1219	void matrix_assign (M &m, const matrix_expression<E> &e) {
1220	typedef R conformant_restrict_type;
1221	typedef typename matrix_assign_traits<typename M::storage_category,
1222	F<typename M::reference, typename E::value_type>::computed,
1223	typename E::const_iterator1::iterator_category,
1224	typename E::const_iterator2::iterator_category>::storage_category storage_category;
1225	// give preference to matrix M's orientation if known
1226	typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>,
1227	typename E::orientation_category ,
1228	typename M::orientation_category >::type orientation_category;
1229	matrix_assign<F, conformant_restrict_type> (m, e, storage_category (), orientation_category ());
1230	}
1231
1232	template<class SC, class RI1, class RI2>
1233	struct matrix_swap_traits {
1234	typedef SC storage_category;
1235	};
1236
1237	template<>
1238	struct matrix_swap_traits<dense_proxy_tag, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
1239	typedef sparse_proxy_tag storage_category;
1240	};
1241
1242	template<>
1243	struct matrix_swap_traits<packed_proxy_tag, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> {
1244	typedef sparse_proxy_tag storage_category;
1245	};
1246
1247	// Dense (proxy) row major case
1248	template<template <class T1, class T2> class F, class R, class M, class E>
1249	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1250	void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, row_major_tag) {
1251	typedef F<typename M::iterator2::reference, typename E::reference> functor_type;
1252	// R unnecessary, make_conformant not required
1253	typedef typename M::size_type size_type;
1254	typedef typename M::difference_type difference_type;
1255	typename M::iterator1 it1 (m.begin1 ());
1256	typename E::iterator1 it1e (e ().begin1 ());
1257	difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (e ().end1 () - it1e)));
1258	while (-- size1 >= 0) {
1259	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1260	typename M::iterator2 it2 (it1.begin ());
1261	typename E::iterator2 it2e (it1e.begin ());
1262	difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (it1e.end () - it2e)));
1263	#else
1264	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1265	typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
1266	difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (end (it1e, iterator1_tag ()) - it2e)));
1267	#endif
1268	while (-- size2 >= 0)
1269	functor_type::apply (it2, it2e), ++ it2, ++ it2e;
1270	++ it1, ++ it1e;
1271	}
1272	}
1273	// Dense (proxy) column major case
1274	template<template <class T1, class T2> class F, class R, class M, class E>
1275	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1276	void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, column_major_tag) {
1277	typedef F<typename M::iterator1::reference, typename E::reference> functor_type;
1278	// R unnecessary, make_conformant not required
1279	typedef typename M::size_type size_type;
1280	typedef typename M::difference_type difference_type;
1281	typename M::iterator2 it2 (m.begin2 ());
1282	typename E::iterator2 it2e (e ().begin2 ());
1283	difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (e ().end2 () - it2e)));
1284	while (-- size2 >= 0) {
1285	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1286	typename M::iterator1 it1 (it2.begin ());
1287	typename E::iterator1 it1e (it2e.begin ());
1288	difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (it2e.end () - it1e)));
1289	#else
1290	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1291	typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
1292	difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (end (it2e, iterator2_tag ()) - it1e)));
1293	#endif
1294	while (-- size1 >= 0)
1295	functor_type::apply (it1, it1e), ++ it1, ++ it1e;
1296	++ it2, ++ it2e;
1297	}
1298	}
1299	// Packed (proxy) row major case
1300	template<template <class T1, class T2> class F, class R, class M, class E>
1301	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1302	void matrix_swap (M &m, matrix_expression<E> &e, packed_proxy_tag, row_major_tag) {
1303	typedef F<typename M::iterator2::reference, typename E::reference> functor_type;
1304	// R unnecessary, make_conformant not required
1305	typedef typename M::size_type size_type;
1306	typedef typename M::difference_type difference_type;
1307	typename M::iterator1 it1 (m.begin1 ());
1308	typename E::iterator1 it1e (e ().begin1 ());
1309	difference_type size1 (BOOST_UBLAS_SAME (m.end1 () - it1, e ().end1 () - it1e));
1310	while (-- size1 >= 0) {
1311	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1312	typename M::iterator2 it2 (it1.begin ());
1313	typename E::iterator2 it2e (it1e.begin ());
1314	difference_type size2 (BOOST_UBLAS_SAME (it1.end () - it2, it1e.end () - it2e));
1315	#else
1316	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1317	typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
1318	difference_type size2 (BOOST_UBLAS_SAME (end (it1, iterator1_tag ()) - it2, end (it1e, iterator1_tag ()) - it2e));
1319	#endif
1320	while (-- size2 >= 0)
1321	functor_type::apply (it2, it2e), ++ it2, ++ it2e;
1322	++ it1, ++ it1e;
1323	}
1324	}
1325	// Packed (proxy) column major case
1326	template<template <class T1, class T2> class F, class R, class M, class E>
1327	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1328	void matrix_swap (M &m, matrix_expression<E> &e, packed_proxy_tag, column_major_tag) {
1329	typedef F<typename M::iterator1::reference, typename E::reference> functor_type;
1330	// R unnecessary, make_conformant not required
1331	typedef typename M::size_type size_type;
1332	typedef typename M::difference_type difference_type;
1333	typename M::iterator2 it2 (m.begin2 ());
1334	typename E::iterator2 it2e (e ().begin2 ());
1335	difference_type size2 (BOOST_UBLAS_SAME (m.end2 () - it2, e ().end2 () - it2e));
1336	while (-- size2 >= 0) {
1337	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1338	typename M::iterator1 it1 (it2.begin ());
1339	typename E::iterator1 it1e (it2e.begin ());
1340	difference_type size1 (BOOST_UBLAS_SAME (it2.end () - it1, it2e.end () - it1e));
1341	#else
1342	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1343	typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
1344	difference_type size1 (BOOST_UBLAS_SAME (end (it2, iterator2_tag ()) - it1, end (it2e, iterator2_tag ()) - it1e));
1345	#endif
1346	while (-- size1 >= 0)
1347	functor_type::apply (it1, it1e), ++ it1, ++ it1e;
1348	++ it2, ++ it2e;
1349	}
1350	}
1351	// Sparse (proxy) row major case
1352	template<template <class T1, class T2> class F, class R, class M, class E>
1353	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1354	void matrix_swap (M &m, matrix_expression<E> &e, sparse_proxy_tag, row_major_tag) {
1355	typedef F<typename M::iterator2::reference, typename E::reference> functor_type;
1356	typedef R conformant_restrict_type;
1357	typedef typename M::size_type size_type;
1358	typedef typename M::difference_type difference_type;
1359	typedef typename M::value_type value_type;
1360	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
1361	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
1362
1363	detail::make_conformant (m, e, row_major_tag (), conformant_restrict_type ());
1364	// FIXME should be a seperate restriction for E
1365	detail::make_conformant (e (), m, row_major_tag (), conformant_restrict_type ());
1366
1367	typename M::iterator1 it1 (m.begin1 ());
1368	typename M::iterator1 it1_end (m.end1 ());
1369	typename E::iterator1 it1e (e ().begin1 ());
1370	typename E::iterator1 it1e_end (e ().end1 ());
1371	while (it1 != it1_end && it1e != it1e_end) {
1372	difference_type compare = it1.index1 () - it1e.index1 ();
1373	if (compare == 0) {
1374	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1375	typename M::iterator2 it2 (it1.begin ());
1376	typename M::iterator2 it2_end (it1.end ());
1377	typename E::iterator2 it2e (it1e.begin ());
1378	typename E::iterator2 it2e_end (it1e.end ());
1379	#else
1380	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1381	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
1382	typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
1383	typename E::iterator2 it2e_end (end (it1e, iterator1_tag ()));
1384	#endif
1385	if (it2 != it2_end && it2e != it2e_end) {
1386	size_type it2_index = it2.index2 (), it2e_index = it2e.index2 ();
1387	while (true) {
1388	difference_type compare = it2_index - it2e_index;
1389	if (compare == 0) {
1390	functor_type::apply (it2, it2e);
1391	++ it2, ++ it2e;
1392	if (it2 != it2_end && it2e != it2e_end) {
1393	it2_index = it2.index2 ();
1394	it2e_index = it2e.index2 ();
1395	} else
1396	break;
1397	} else if (compare < 0) {
1398	increment (it2, it2_end, - compare);
1399	if (it2 != it2_end)
1400	it2_index = it2.index2 ();
1401	else
1402	break;
1403	} else if (compare > 0) {
1404	increment (it2e, it2e_end, compare);
1405	if (it2e != it2e_end)
1406	it2e_index = it2e.index2 ();
1407	else
1408	break;
1409	}
1410	}
1411	}
1412	#if BOOST_UBLAS_TYPE_CHECK
1413	increment (it2e, it2e_end);
1414	increment (it2, it2_end);
1415	#endif
1416	++ it1, ++ it1e;
1417	} else if (compare < 0) {
1418	#if BOOST_UBLAS_TYPE_CHECK
1419	while (it1.index1 () < it1e.index1 ()) {
1420	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1421	typename M::iterator2 it2 (it1.begin ());
1422	typename M::iterator2 it2_end (it1.end ());
1423	#else
1424	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1425	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
1426	#endif
1427	increment (it2, it2_end);
1428	++ it1;
1429	}
1430	#else
1431	increment (it1, it1_end, - compare);
1432	#endif
1433	} else if (compare > 0) {
1434	#if BOOST_UBLAS_TYPE_CHECK
1435	while (it1e.index1 () < it1.index1 ()) {
1436	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1437	typename E::iterator2 it2e (it1e.begin ());
1438	typename E::iterator2 it2e_end (it1e.end ());
1439	#else
1440	typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
1441	typename E::iterator2 it2e_end (end (it1e, iterator1_tag ()));
1442	#endif
1443	increment (it2e, it2e_end);
1444	++ it1e;
1445	}
1446	#else
1447	increment (it1e, it1e_end, compare);
1448	#endif
1449	}
1450	}
1451	#if BOOST_UBLAS_TYPE_CHECK
1452	while (it1e != it1e_end) {
1453	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1454	typename E::iterator2 it2e (it1e.begin ());
1455	typename E::iterator2 it2e_end (it1e.end ());
1456	#else
1457	typename E::iterator2 it2e (begin (it1e, iterator1_tag ()));
1458	typename E::iterator2 it2e_end (end (it1e, iterator1_tag ()));
1459	#endif
1460	increment (it2e, it2e_end);
1461	++ it1e;
1462	}
1463	while (it1 != it1_end) {
1464	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1465	typename M::iterator2 it2 (it1.begin ());
1466	typename M::iterator2 it2_end (it1.end ());
1467	#else
1468	typename M::iterator2 it2 (begin (it1, iterator1_tag ()));
1469	typename M::iterator2 it2_end (end (it1, iterator1_tag ()));
1470	#endif
1471	increment (it2, it2_end);
1472	++ it1;
1473	}
1474	#endif
1475	}
1476	// Sparse (proxy) column major case
1477	template<template <class T1, class T2> class F, class R, class M, class E>
1478	// BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it.
1479	void matrix_swap (M &m, matrix_expression<E> &e, sparse_proxy_tag, column_major_tag) {
1480	typedef F<typename M::iterator1::reference, typename E::reference> functor_type;
1481	typedef R conformant_restrict_type;
1482	typedef typename M::size_type size_type;
1483	typedef typename M::difference_type difference_type;
1484	typedef typename M::value_type value_type;
1485	BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ());
1486	BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ());
1487
1488	detail::make_conformant (m, e, column_major_tag (), conformant_restrict_type ());
1489	// FIXME should be a seperate restriction for E
1490	detail::make_conformant (e (), m, column_major_tag (), conformant_restrict_type ());
1491
1492	typename M::iterator2 it2 (m.begin2 ());
1493	typename M::iterator2 it2_end (m.end2 ());
1494	typename E::iterator2 it2e (e ().begin2 ());
1495	typename E::iterator2 it2e_end (e ().end2 ());
1496	while (it2 != it2_end && it2e != it2e_end) {
1497	difference_type compare = it2.index2 () - it2e.index2 ();
1498	if (compare == 0) {
1499	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1500	typename M::iterator1 it1 (it2.begin ());
1501	typename M::iterator1 it1_end (it2.end ());
1502	typename E::iterator1 it1e (it2e.begin ());
1503	typename E::iterator1 it1e_end (it2e.end ());
1504	#else
1505	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1506	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1507	typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
1508	typename E::iterator1 it1e_end (end (it2e, iterator2_tag ()));
1509	#endif
1510	if (it1 != it1_end && it1e != it1e_end) {
1511	size_type it1_index = it1.index1 (), it1e_index = it1e.index1 ();
1512	while (true) {
1513	difference_type compare = it1_index - it1e_index;
1514	if (compare == 0) {
1515	functor_type::apply (it1, it1e);
1516	++ it1, ++ it1e;
1517	if (it1 != it1_end && it1e != it1e_end) {
1518	it1_index = it1.index1 ();
1519	it1e_index = it1e.index1 ();
1520	} else
1521	break;
1522	} else if (compare < 0) {
1523	increment (it1, it1_end, - compare);
1524	if (it1 != it1_end)
1525	it1_index = it1.index1 ();
1526	else
1527	break;
1528	} else if (compare > 0) {
1529	increment (it1e, it1e_end, compare);
1530	if (it1e != it1e_end)
1531	it1e_index = it1e.index1 ();
1532	else
1533	break;
1534	}
1535	}
1536	}
1537	#if BOOST_UBLAS_TYPE_CHECK
1538	increment (it1e, it1e_end);
1539	increment (it1, it1_end);
1540	#endif
1541	++ it2, ++ it2e;
1542	} else if (compare < 0) {
1543	#if BOOST_UBLAS_TYPE_CHECK
1544	while (it2.index2 () < it2e.index2 ()) {
1545	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1546	typename M::iterator1 it1 (it2.begin ());
1547	typename M::iterator1 it1_end (it2.end ());
1548	#else
1549	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1550	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1551	#endif
1552	increment (it1, it1_end);
1553	++ it2;
1554	}
1555	#else
1556	increment (it2, it2_end, - compare);
1557	#endif
1558	} else if (compare > 0) {
1559	#if BOOST_UBLAS_TYPE_CHECK
1560	while (it2e.index2 () < it2.index2 ()) {
1561	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1562	typename E::iterator1 it1e (it2e.begin ());
1563	typename E::iterator1 it1e_end (it2e.end ());
1564	#else
1565	typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
1566	typename E::iterator1 it1e_end (end (it2e, iterator2_tag ()));
1567	#endif
1568	increment (it1e, it1e_end);
1569	++ it2e;
1570	}
1571	#else
1572	increment (it2e, it2e_end, compare);
1573	#endif
1574	}
1575	}
1576	#if BOOST_UBLAS_TYPE_CHECK
1577	while (it2e != it2e_end) {
1578	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1579	typename E::iterator1 it1e (it2e.begin ());
1580	typename E::iterator1 it1e_end (it2e.end ());
1581	#else
1582	typename E::iterator1 it1e (begin (it2e, iterator2_tag ()));
1583	typename E::iterator1 it1e_end (end (it2e, iterator2_tag ()));
1584	#endif
1585	increment (it1e, it1e_end);
1586	++ it2e;
1587	}
1588	while (it2 != it2_end) {
1589	#ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION
1590	typename M::iterator1 it1 (it2.begin ());
1591	typename M::iterator1 it1_end (it2.end ());
1592	#else
1593	typename M::iterator1 it1 (begin (it2, iterator2_tag ()));
1594	typename M::iterator1 it1_end (end (it2, iterator2_tag ()));
1595	#endif
1596	increment (it1, it1_end);
1597	++ it2;
1598	}
1599	#endif
1600	}
1601
1602	// Dispatcher
1603	template<template <class T1, class T2> class F, class M, class E>
1604	BOOST_UBLAS_INLINE
1605	void matrix_swap (M &m, matrix_expression<E> &e) {
1606	typedef typename matrix_swap_traits<typename M::storage_category,
1607	typename E::const_iterator1::iterator_category,
1608	typename E::const_iterator2::iterator_category>::storage_category storage_category;
1609	// give preference to matrix M's orientation if known
1610	typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>,
1611	typename E::orientation_category ,
1612	typename M::orientation_category >::type orientation_category;
1613	typedef basic_full<typename M::size_type> unrestricted;
1614	matrix_swap<F, unrestricted> (m, e, storage_category (), orientation_category ());
1615	}
1616	template<template <class T1, class T2> class F, class R, class M, class E>
1617	BOOST_UBLAS_INLINE
1618	void matrix_swap (M &m, matrix_expression<E> &e) {
1619	typedef R conformant_restrict_type;
1620	typedef typename matrix_swap_traits<typename M::storage_category,
1621	typename E::const_iterator1::iterator_category,
1622	typename E::const_iterator2::iterator_category>::storage_category storage_category;
1623	// give preference to matrix M's orientation if known
1624	typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>,
1625	typename E::orientation_category ,
1626	typename M::orientation_category >::type orientation_category;
1627	matrix_swap<F, conformant_restrict_type> (m, e, storage_category (), orientation_category ());
1628	}
1629
1630	}}}
1631
1632	#endif

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

Context Navigation

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

Download in other formats: