![]() |
Boost.uBlas 1.49
Linear Algebra in C++: matrices, vectors and numeric algorithms
|
00001 // 00002 // Copyright (c) 2000-2002 00003 // Joerg Walter, Mathias Koch 00004 // 00005 // Distributed under the Boost Software License, Version 1.0. (See 00006 // accompanying file LICENSE_1_0.txt or copy at 00007 // http://www.boost.org/LICENSE_1_0.txt) 00008 // 00009 // The authors gratefully acknowledge the support of 00010 // GeNeSys mbH & Co. KG in producing this work. 00011 // 00012 00013 #ifndef _BOOST_UBLAS_MATRIX_ASSIGN_ 00014 #define _BOOST_UBLAS_MATRIX_ASSIGN_ 00015 00016 #include <boost/numeric/ublas/traits.hpp> 00017 // Required for make_conformant storage 00018 #include <vector> 00019 00020 // Iterators based on ideas of Jeremy Siek 00021 00022 namespace boost { namespace numeric { namespace ublas { 00023 namespace detail { 00024 00025 // Weak equality check - useful to compare equality two arbitary matrix expression results. 00026 // Since the actual expressions are unknown, we check for and arbitary error bound 00027 // on the relative error. 00028 // For a linear expression the infinity norm makes sense as we do not know how the elements will be 00029 // combined in the expression. False positive results are inevitable for arbirary expressions! 00030 template<class E1, class E2, class S> 00031 BOOST_UBLAS_INLINE 00032 bool equals (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2, S epsilon, S min_norm) { 00033 return norm_inf (e1 - e2) < epsilon * 00034 std::max<S> (std::max<S> (norm_inf (e1), norm_inf (e2)), min_norm); 00035 } 00036 00037 template<class E1, class E2> 00038 BOOST_UBLAS_INLINE 00039 bool expression_type_check (const matrix_expression<E1> &e1, const matrix_expression<E2> &e2) { 00040 typedef typename type_traits<typename promote_traits<typename E1::value_type, 00041 typename E2::value_type>::promote_type>::real_type real_type; 00042 return equals (e1, e2, BOOST_UBLAS_TYPE_CHECK_EPSILON, BOOST_UBLAS_TYPE_CHECK_MIN); 00043 } 00044 00045 00046 template<class M, class E, class R> 00047 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00048 void make_conformant (M &m, const matrix_expression<E> &e, row_major_tag, R) { 00049 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00050 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00051 typedef R conformant_restrict_type; 00052 typedef typename M::size_type size_type; 00053 typedef typename M::difference_type difference_type; 00054 typedef typename M::value_type value_type; 00055 // FIXME unbounded_array with push_back maybe better 00056 std::vector<std::pair<size_type, size_type> > index; 00057 typename M::iterator1 it1 (m.begin1 ()); 00058 typename M::iterator1 it1_end (m.end1 ()); 00059 typename E::const_iterator1 it1e (e ().begin1 ()); 00060 typename E::const_iterator1 it1e_end (e ().end1 ()); 00061 while (it1 != it1_end && it1e != it1e_end) { 00062 difference_type compare = it1.index1 () - it1e.index1 (); 00063 if (compare == 0) { 00064 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00065 typename M::iterator2 it2 (it1.begin ()); 00066 typename M::iterator2 it2_end (it1.end ()); 00067 typename E::const_iterator2 it2e (it1e.begin ()); 00068 typename E::const_iterator2 it2e_end (it1e.end ()); 00069 #else 00070 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00071 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00072 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00073 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00074 #endif 00075 if (it2 != it2_end && it2e != it2e_end) { 00076 size_type it2_index = it2.index2 (), it2e_index = it2e.index2 (); 00077 while (true) { 00078 difference_type compare = it2_index - it2e_index; 00079 if (compare == 0) { 00080 ++ it2, ++ it2e; 00081 if (it2 != it2_end && it2e != it2e_end) { 00082 it2_index = it2.index2 (); 00083 it2e_index = it2e.index2 (); 00084 } else 00085 break; 00086 } else if (compare < 0) { 00087 increment (it2, it2_end, - compare); 00088 if (it2 != it2_end) 00089 it2_index = it2.index2 (); 00090 else 00091 break; 00092 } else if (compare > 0) { 00093 if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ())) 00094 if (static_cast<value_type>(*it2e) != value_type/*zero*/()) 00095 index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ())); 00096 ++ it2e; 00097 if (it2e != it2e_end) 00098 it2e_index = it2e.index2 (); 00099 else 00100 break; 00101 } 00102 } 00103 } 00104 while (it2e != it2e_end) { 00105 if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ())) 00106 if (static_cast<value_type>(*it2e) != value_type/*zero*/()) 00107 index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ())); 00108 ++ it2e; 00109 } 00110 ++ it1, ++ it1e; 00111 } else if (compare < 0) { 00112 increment (it1, it1_end, - compare); 00113 } else if (compare > 0) { 00114 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00115 typename E::const_iterator2 it2e (it1e.begin ()); 00116 typename E::const_iterator2 it2e_end (it1e.end ()); 00117 #else 00118 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00119 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00120 #endif 00121 while (it2e != it2e_end) { 00122 if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ())) 00123 if (static_cast<value_type>(*it2e) != value_type/*zero*/()) 00124 index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ())); 00125 ++ it2e; 00126 } 00127 ++ it1e; 00128 } 00129 } 00130 while (it1e != it1e_end) { 00131 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00132 typename E::const_iterator2 it2e (it1e.begin ()); 00133 typename E::const_iterator2 it2e_end (it1e.end ()); 00134 #else 00135 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00136 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00137 #endif 00138 while (it2e != it2e_end) { 00139 if (conformant_restrict_type::other (it2e.index1 (), it2e.index2 ())) 00140 if (static_cast<value_type>(*it2e) != value_type/*zero*/()) 00141 index.push_back (std::pair<size_type, size_type> (it2e.index1 (), it2e.index2 ())); 00142 ++ it2e; 00143 } 00144 ++ it1e; 00145 } 00146 // ISSUE proxies require insert_element 00147 for (size_type k = 0; k < index.size (); ++ k) 00148 m (index [k].first, index [k].second) = value_type/*zero*/(); 00149 } 00150 template<class M, class E, class R> 00151 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00152 void make_conformant (M &m, const matrix_expression<E> &e, column_major_tag, R) { 00153 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00154 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00155 typedef R conformant_restrict_type; 00156 typedef typename M::size_type size_type; 00157 typedef typename M::difference_type difference_type; 00158 typedef typename M::value_type value_type; 00159 std::vector<std::pair<size_type, size_type> > index; 00160 typename M::iterator2 it2 (m.begin2 ()); 00161 typename M::iterator2 it2_end (m.end2 ()); 00162 typename E::const_iterator2 it2e (e ().begin2 ()); 00163 typename E::const_iterator2 it2e_end (e ().end2 ()); 00164 while (it2 != it2_end && it2e != it2e_end) { 00165 difference_type compare = it2.index2 () - it2e.index2 (); 00166 if (compare == 0) { 00167 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00168 typename M::iterator1 it1 (it2.begin ()); 00169 typename M::iterator1 it1_end (it2.end ()); 00170 typename E::const_iterator1 it1e (it2e.begin ()); 00171 typename E::const_iterator1 it1e_end (it2e.end ()); 00172 #else 00173 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00174 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 00175 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00176 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 00177 #endif 00178 if (it1 != it1_end && it1e != it1e_end) { 00179 size_type it1_index = it1.index1 (), it1e_index = it1e.index1 (); 00180 while (true) { 00181 difference_type compare = it1_index - it1e_index; 00182 if (compare == 0) { 00183 ++ it1, ++ it1e; 00184 if (it1 != it1_end && it1e != it1e_end) { 00185 it1_index = it1.index1 (); 00186 it1e_index = it1e.index1 (); 00187 } else 00188 break; 00189 } else if (compare < 0) { 00190 increment (it1, it1_end, - compare); 00191 if (it1 != it1_end) 00192 it1_index = it1.index1 (); 00193 else 00194 break; 00195 } else if (compare > 0) { 00196 if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ())) 00197 if (static_cast<value_type>(*it1e) != value_type/*zero*/()) 00198 index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ())); 00199 ++ it1e; 00200 if (it1e != it1e_end) 00201 it1e_index = it1e.index1 (); 00202 else 00203 break; 00204 } 00205 } 00206 } 00207 while (it1e != it1e_end) { 00208 if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ())) 00209 if (static_cast<value_type>(*it1e) != value_type/*zero*/()) 00210 index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ())); 00211 ++ it1e; 00212 } 00213 ++ it2, ++ it2e; 00214 } else if (compare < 0) { 00215 increment (it2, it2_end, - compare); 00216 } else if (compare > 0) { 00217 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00218 typename E::const_iterator1 it1e (it2e.begin ()); 00219 typename E::const_iterator1 it1e_end (it2e.end ()); 00220 #else 00221 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00222 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 00223 #endif 00224 while (it1e != it1e_end) { 00225 if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ())) 00226 if (static_cast<value_type>(*it1e) != value_type/*zero*/()) 00227 index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ())); 00228 ++ it1e; 00229 } 00230 ++ it2e; 00231 } 00232 } 00233 while (it2e != it2e_end) { 00234 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00235 typename E::const_iterator1 it1e (it2e.begin ()); 00236 typename E::const_iterator1 it1e_end (it2e.end ()); 00237 #else 00238 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00239 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 00240 #endif 00241 while (it1e != it1e_end) { 00242 if (conformant_restrict_type::other (it1e.index1 (), it1e.index2 ())) 00243 if (static_cast<value_type>(*it1e) != value_type/*zero*/()) 00244 index.push_back (std::pair<size_type, size_type> (it1e.index1 (), it1e.index2 ())); 00245 ++ it1e; 00246 } 00247 ++ it2e; 00248 } 00249 // ISSUE proxies require insert_element 00250 for (size_type k = 0; k < index.size (); ++ k) 00251 m (index [k].first, index [k].second) = value_type/*zero*/(); 00252 } 00253 00254 }//namespace detail 00255 00256 00257 // Explicitly iterating row major 00258 template<template <class T1, class T2> class F, class M, class T> 00259 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00260 void iterating_matrix_assign_scalar (M &m, const T &t, row_major_tag) { 00261 typedef F<typename M::iterator2::reference, T> functor_type; 00262 typedef typename M::difference_type difference_type; 00263 difference_type size1 (m.size1 ()); 00264 difference_type size2 (m.size2 ()); 00265 typename M::iterator1 it1 (m.begin1 ()); 00266 BOOST_UBLAS_CHECK (size2 == 0 || m.end1 () - it1 == size1, bad_size ()); 00267 while (-- size1 >= 0) { 00268 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00269 typename M::iterator2 it2 (it1.begin ()); 00270 #else 00271 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00272 #endif 00273 BOOST_UBLAS_CHECK (it1.end () - it2 == size2, bad_size ()); 00274 difference_type temp_size2 (size2); 00275 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00276 while (-- temp_size2 >= 0) 00277 functor_type::apply (*it2, t), ++ it2; 00278 #else 00279 DD (temp_size2, 4, r, (functor_type::apply (*it2, t), ++ it2)); 00280 #endif 00281 ++ it1; 00282 } 00283 } 00284 // Explicitly iterating column major 00285 template<template <class T1, class T2> class F, class M, class T> 00286 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00287 void iterating_matrix_assign_scalar (M &m, const T &t, column_major_tag) { 00288 typedef F<typename M::iterator1::reference, T> functor_type; 00289 typedef typename M::difference_type difference_type; 00290 difference_type size2 (m.size2 ()); 00291 difference_type size1 (m.size1 ()); 00292 typename M::iterator2 it2 (m.begin2 ()); 00293 BOOST_UBLAS_CHECK (size1 == 0 || m.end2 () - it2 == size2, bad_size ()); 00294 while (-- size2 >= 0) { 00295 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00296 typename M::iterator1 it1 (it2.begin ()); 00297 #else 00298 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00299 #endif 00300 BOOST_UBLAS_CHECK (it2.end () - it1 == size1, bad_size ()); 00301 difference_type temp_size1 (size1); 00302 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00303 while (-- temp_size1 >= 0) 00304 functor_type::apply (*it1, t), ++ it1; 00305 #else 00306 DD (temp_size1, 4, r, (functor_type::apply (*it1, t), ++ it1)); 00307 #endif 00308 ++ it2; 00309 } 00310 } 00311 // Explicitly indexing row major 00312 template<template <class T1, class T2> class F, class M, class T> 00313 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00314 void indexing_matrix_assign_scalar (M &m, const T &t, row_major_tag) { 00315 typedef F<typename M::reference, T> functor_type; 00316 typedef typename M::size_type size_type; 00317 size_type size1 (m.size1 ()); 00318 size_type size2 (m.size2 ()); 00319 for (size_type i = 0; i < size1; ++ i) { 00320 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00321 for (size_type j = 0; j < size2; ++ j) 00322 functor_type::apply (m (i, j), t); 00323 #else 00324 size_type j (0); 00325 DD (size2, 4, r, (functor_type::apply (m (i, j), t), ++ j)); 00326 #endif 00327 } 00328 } 00329 // Explicitly indexing column major 00330 template<template <class T1, class T2> class F, class M, class T> 00331 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00332 void indexing_matrix_assign_scalar (M &m, const T &t, column_major_tag) { 00333 typedef F<typename M::reference, T> functor_type; 00334 typedef typename M::size_type size_type; 00335 size_type size2 (m.size2 ()); 00336 size_type size1 (m.size1 ()); 00337 for (size_type j = 0; j < size2; ++ j) { 00338 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00339 for (size_type i = 0; i < size1; ++ i) 00340 functor_type::apply (m (i, j), t); 00341 #else 00342 size_type i (0); 00343 DD (size1, 4, r, (functor_type::apply (m (i, j), t), ++ i)); 00344 #endif 00345 } 00346 } 00347 00348 // Dense (proxy) case 00349 template<template <class T1, class T2> class F, class M, class T, class C> 00350 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00351 void matrix_assign_scalar (M &m, const T &t, dense_proxy_tag, C) { 00352 typedef C orientation_category; 00353 #ifdef BOOST_UBLAS_USE_INDEXING 00354 indexing_matrix_assign_scalar<F> (m, t, orientation_category ()); 00355 #elif BOOST_UBLAS_USE_ITERATING 00356 iterating_matrix_assign_scalar<F> (m, t, orientation_category ()); 00357 #else 00358 typedef typename M::size_type size_type; 00359 size_type size1 (m.size1 ()); 00360 size_type size2 (m.size2 ()); 00361 if (size1 >= BOOST_UBLAS_ITERATOR_THRESHOLD && 00362 size2 >= BOOST_UBLAS_ITERATOR_THRESHOLD) 00363 iterating_matrix_assign_scalar<F> (m, t, orientation_category ()); 00364 else 00365 indexing_matrix_assign_scalar<F> (m, t, orientation_category ()); 00366 #endif 00367 } 00368 // Packed (proxy) row major case 00369 template<template <class T1, class T2> class F, class M, class T> 00370 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00371 void matrix_assign_scalar (M &m, const T &t, packed_proxy_tag, row_major_tag) { 00372 typedef F<typename M::iterator2::reference, T> functor_type; 00373 typedef typename M::difference_type difference_type; 00374 typename M::iterator1 it1 (m.begin1 ()); 00375 difference_type size1 (m.end1 () - it1); 00376 while (-- size1 >= 0) { 00377 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00378 typename M::iterator2 it2 (it1.begin ()); 00379 difference_type size2 (it1.end () - it2); 00380 #else 00381 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00382 difference_type size2 (end (it1, iterator1_tag ()) - it2); 00383 #endif 00384 while (-- size2 >= 0) 00385 functor_type::apply (*it2, t), ++ it2; 00386 ++ it1; 00387 } 00388 } 00389 // Packed (proxy) column major case 00390 template<template <class T1, class T2> class F, class M, class T> 00391 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00392 void matrix_assign_scalar (M &m, const T &t, packed_proxy_tag, column_major_tag) { 00393 typedef F<typename M::iterator1::reference, T> functor_type; 00394 typedef typename M::difference_type difference_type; 00395 typename M::iterator2 it2 (m.begin2 ()); 00396 difference_type size2 (m.end2 () - it2); 00397 while (-- size2 >= 0) { 00398 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00399 typename M::iterator1 it1 (it2.begin ()); 00400 difference_type size1 (it2.end () - it1); 00401 #else 00402 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00403 difference_type size1 (end (it2, iterator2_tag ()) - it1); 00404 #endif 00405 while (-- size1 >= 0) 00406 functor_type::apply (*it1, t), ++ it1; 00407 ++ it2; 00408 } 00409 } 00410 // Sparse (proxy) row major case 00411 template<template <class T1, class T2> class F, class M, class T> 00412 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00413 void matrix_assign_scalar (M &m, const T &t, sparse_proxy_tag, row_major_tag) { 00414 typedef F<typename M::iterator2::reference, T> functor_type; 00415 typename M::iterator1 it1 (m.begin1 ()); 00416 typename M::iterator1 it1_end (m.end1 ()); 00417 while (it1 != it1_end) { 00418 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00419 typename M::iterator2 it2 (it1.begin ()); 00420 typename M::iterator2 it2_end (it1.end ()); 00421 #else 00422 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00423 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00424 #endif 00425 while (it2 != it2_end) 00426 functor_type::apply (*it2, t), ++ it2; 00427 ++ it1; 00428 } 00429 } 00430 // Sparse (proxy) column major case 00431 template<template <class T1, class T2> class F, class M, class T> 00432 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00433 void matrix_assign_scalar (M &m, const T &t, sparse_proxy_tag, column_major_tag) { 00434 typedef F<typename M::iterator1::reference, T> functor_type; 00435 typename M::iterator2 it2 (m.begin2 ()); 00436 typename M::iterator2 it2_end (m.end2 ()); 00437 while (it2 != it2_end) { 00438 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00439 typename M::iterator1 it1 (it2.begin ()); 00440 typename M::iterator1 it1_end (it2.end ()); 00441 #else 00442 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00443 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 00444 #endif 00445 while (it1 != it1_end) 00446 functor_type::apply (*it1, t), ++ it1; 00447 ++ it2; 00448 } 00449 } 00450 00451 // Dispatcher 00452 template<template <class T1, class T2> class F, class M, class T> 00453 BOOST_UBLAS_INLINE 00454 void matrix_assign_scalar (M &m, const T &t) { 00455 typedef typename M::storage_category storage_category; 00456 typedef typename M::orientation_category orientation_category; 00457 matrix_assign_scalar<F> (m, t, storage_category (), orientation_category ()); 00458 } 00459 00460 template<class SC, bool COMPUTED, class RI1, class RI2> 00461 struct matrix_assign_traits { 00462 typedef SC storage_category; 00463 }; 00464 00465 template<bool COMPUTED> 00466 struct matrix_assign_traits<dense_tag, COMPUTED, packed_random_access_iterator_tag, packed_random_access_iterator_tag> { 00467 typedef packed_tag storage_category; 00468 }; 00469 template<> 00470 struct matrix_assign_traits<dense_tag, false, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00471 typedef sparse_tag storage_category; 00472 }; 00473 template<> 00474 struct matrix_assign_traits<dense_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00475 typedef sparse_proxy_tag storage_category; 00476 }; 00477 00478 template<bool COMPUTED> 00479 struct matrix_assign_traits<dense_proxy_tag, COMPUTED, packed_random_access_iterator_tag, packed_random_access_iterator_tag> { 00480 typedef packed_proxy_tag storage_category; 00481 }; 00482 template<bool COMPUTED> 00483 struct matrix_assign_traits<dense_proxy_tag, COMPUTED, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00484 typedef sparse_proxy_tag storage_category; 00485 }; 00486 00487 template<> 00488 struct matrix_assign_traits<packed_tag, false, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00489 typedef sparse_tag storage_category; 00490 }; 00491 template<> 00492 struct matrix_assign_traits<packed_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00493 typedef sparse_proxy_tag storage_category; 00494 }; 00495 00496 template<bool COMPUTED> 00497 struct matrix_assign_traits<packed_proxy_tag, COMPUTED, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00498 typedef sparse_proxy_tag storage_category; 00499 }; 00500 00501 template<> 00502 struct matrix_assign_traits<sparse_tag, true, dense_random_access_iterator_tag, dense_random_access_iterator_tag> { 00503 typedef sparse_proxy_tag storage_category; 00504 }; 00505 template<> 00506 struct matrix_assign_traits<sparse_tag, true, packed_random_access_iterator_tag, packed_random_access_iterator_tag> { 00507 typedef sparse_proxy_tag storage_category; 00508 }; 00509 template<> 00510 struct matrix_assign_traits<sparse_tag, true, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 00511 typedef sparse_proxy_tag storage_category; 00512 }; 00513 00514 // Explicitly iterating row major 00515 template<template <class T1, class T2> class F, class M, class E> 00516 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00517 void iterating_matrix_assign (M &m, const matrix_expression<E> &e, row_major_tag) { 00518 typedef F<typename M::iterator2::reference, typename E::value_type> functor_type; 00519 typedef typename M::difference_type difference_type; 00520 difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ())); 00521 difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ())); 00522 typename M::iterator1 it1 (m.begin1 ()); 00523 BOOST_UBLAS_CHECK (size2 == 0 || m.end1 () - it1 == size1, bad_size ()); 00524 typename E::const_iterator1 it1e (e ().begin1 ()); 00525 BOOST_UBLAS_CHECK (size2 == 0 || e ().end1 () - it1e == size1, bad_size ()); 00526 while (-- size1 >= 0) { 00527 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00528 typename M::iterator2 it2 (it1.begin ()); 00529 typename E::const_iterator2 it2e (it1e.begin ()); 00530 #else 00531 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00532 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00533 #endif 00534 BOOST_UBLAS_CHECK (it1.end () - it2 == size2, bad_size ()); 00535 BOOST_UBLAS_CHECK (it1e.end () - it2e == size2, bad_size ()); 00536 difference_type temp_size2 (size2); 00537 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00538 while (-- temp_size2 >= 0) 00539 functor_type::apply (*it2, *it2e), ++ it2, ++ it2e; 00540 #else 00541 DD (temp_size2, 2, r, (functor_type::apply (*it2, *it2e), ++ it2, ++ it2e)); 00542 #endif 00543 ++ it1, ++ it1e; 00544 } 00545 } 00546 // Explicitly iterating column major 00547 template<template <class T1, class T2> class F, class M, class E> 00548 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00549 void iterating_matrix_assign (M &m, const matrix_expression<E> &e, column_major_tag) { 00550 typedef F<typename M::iterator1::reference, typename E::value_type> functor_type; 00551 typedef typename M::difference_type difference_type; 00552 difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ())); 00553 difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ())); 00554 typename M::iterator2 it2 (m.begin2 ()); 00555 BOOST_UBLAS_CHECK (size1 == 0 || m.end2 () - it2 == size2, bad_size ()); 00556 typename E::const_iterator2 it2e (e ().begin2 ()); 00557 BOOST_UBLAS_CHECK (size1 == 0 || e ().end2 () - it2e == size2, bad_size ()); 00558 while (-- size2 >= 0) { 00559 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00560 typename M::iterator1 it1 (it2.begin ()); 00561 typename E::const_iterator1 it1e (it2e.begin ()); 00562 #else 00563 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00564 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00565 #endif 00566 BOOST_UBLAS_CHECK (it2.end () - it1 == size1, bad_size ()); 00567 BOOST_UBLAS_CHECK (it2e.end () - it1e == size1, bad_size ()); 00568 difference_type temp_size1 (size1); 00569 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00570 while (-- temp_size1 >= 0) 00571 functor_type::apply (*it1, *it1e), ++ it1, ++ it1e; 00572 #else 00573 DD (temp_size1, 2, r, (functor_type::apply (*it1, *it1e), ++ it1, ++ it1e)); 00574 #endif 00575 ++ it2, ++ it2e; 00576 } 00577 } 00578 // Explicitly indexing row major 00579 template<template <class T1, class T2> class F, class M, class E> 00580 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00581 void indexing_matrix_assign (M &m, const matrix_expression<E> &e, row_major_tag) { 00582 typedef F<typename M::reference, typename E::value_type> functor_type; 00583 typedef typename M::size_type size_type; 00584 size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ())); 00585 size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ())); 00586 for (size_type i = 0; i < size1; ++ i) { 00587 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00588 for (size_type j = 0; j < size2; ++ j) 00589 functor_type::apply (m (i, j), e () (i, j)); 00590 #else 00591 size_type j (0); 00592 DD (size2, 2, r, (functor_type::apply (m (i, j), e () (i, j)), ++ j)); 00593 #endif 00594 } 00595 } 00596 // Explicitly indexing column major 00597 template<template <class T1, class T2> class F, class M, class E> 00598 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00599 void indexing_matrix_assign (M &m, const matrix_expression<E> &e, column_major_tag) { 00600 typedef F<typename M::reference, typename E::value_type> functor_type; 00601 typedef typename M::size_type size_type; 00602 size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ())); 00603 size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ())); 00604 for (size_type j = 0; j < size2; ++ j) { 00605 #ifndef BOOST_UBLAS_USE_DUFF_DEVICE 00606 for (size_type i = 0; i < size1; ++ i) 00607 functor_type::apply (m (i, j), e () (i, j)); 00608 #else 00609 size_type i (0); 00610 DD (size1, 2, r, (functor_type::apply (m (i, j), e () (i, j)), ++ i)); 00611 #endif 00612 } 00613 } 00614 00615 // Dense (proxy) case 00616 template<template <class T1, class T2> class F, class R, class M, class E, class C> 00617 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00618 void matrix_assign (M &m, const matrix_expression<E> &e, dense_proxy_tag, C) { 00619 // R unnecessary, make_conformant not required 00620 typedef C orientation_category; 00621 #ifdef BOOST_UBLAS_USE_INDEXING 00622 indexing_matrix_assign<F> (m, e, orientation_category ()); 00623 #elif BOOST_UBLAS_USE_ITERATING 00624 iterating_matrix_assign<F> (m, e, orientation_category ()); 00625 #else 00626 typedef typename M::difference_type difference_type; 00627 size_type size1 (BOOST_UBLAS_SAME (m.size1 (), e ().size1 ())); 00628 size_type size2 (BOOST_UBLAS_SAME (m.size2 (), e ().size2 ())); 00629 if (size1 >= BOOST_UBLAS_ITERATOR_THRESHOLD && 00630 size2 >= BOOST_UBLAS_ITERATOR_THRESHOLD) 00631 iterating_matrix_assign<F> (m, e, orientation_category ()); 00632 else 00633 indexing_matrix_assign<F> (m, e, orientation_category ()); 00634 #endif 00635 } 00636 // Packed (proxy) row major case 00637 template<template <class T1, class T2> class F, class R, class M, class E> 00638 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00639 void matrix_assign (M &m, const matrix_expression<E> &e, packed_proxy_tag, row_major_tag) { 00640 typedef typename matrix_traits<E>::value_type expr_value_type; 00641 typedef F<typename M::iterator2::reference, expr_value_type> functor_type; 00642 // R unnecessary, make_conformant not required 00643 typedef typename M::difference_type difference_type; 00644 typedef typename M::value_type value_type; 00645 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00646 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00647 #if BOOST_UBLAS_TYPE_CHECK 00648 matrix<value_type, row_major> cm (m.size1 (), m.size2 ()); 00649 indexing_matrix_assign<scalar_assign> (cm, m, row_major_tag ()); 00650 indexing_matrix_assign<F> (cm, e, row_major_tag ()); 00651 #endif 00652 typename M::iterator1 it1 (m.begin1 ()); 00653 typename M::iterator1 it1_end (m.end1 ()); 00654 typename E::const_iterator1 it1e (e ().begin1 ()); 00655 typename E::const_iterator1 it1e_end (e ().end1 ()); 00656 difference_type it1_size (it1_end - it1); 00657 difference_type it1e_size (it1e_end - it1e); 00658 difference_type diff1 (0); 00659 if (it1_size > 0 && it1e_size > 0) 00660 diff1 = it1.index1 () - it1e.index1 (); 00661 if (diff1 != 0) { 00662 difference_type size1 = (std::min) (diff1, it1e_size); 00663 if (size1 > 0) { 00664 it1e += size1; 00665 it1e_size -= size1; 00666 diff1 -= size1; 00667 } 00668 size1 = (std::min) (- diff1, it1_size); 00669 if (size1 > 0) { 00670 it1_size -= size1; 00671 if (!functor_type::computed) { 00672 while (-- size1 >= 0) { // zeroing 00673 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00674 typename M::iterator2 it2 (it1.begin ()); 00675 typename M::iterator2 it2_end (it1.end ()); 00676 #else 00677 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00678 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00679 #endif 00680 difference_type size2 (it2_end - it2); 00681 while (-- size2 >= 0) 00682 functor_type::apply (*it2, expr_value_type/*zero*/()), ++ it2; 00683 ++ it1; 00684 } 00685 } else { 00686 it1 += size1; 00687 } 00688 diff1 += size1; 00689 } 00690 } 00691 difference_type size1 ((std::min) (it1_size, it1e_size)); 00692 it1_size -= size1; 00693 it1e_size -= size1; 00694 while (-- size1 >= 0) { 00695 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00696 typename M::iterator2 it2 (it1.begin ()); 00697 typename M::iterator2 it2_end (it1.end ()); 00698 typename E::const_iterator2 it2e (it1e.begin ()); 00699 typename E::const_iterator2 it2e_end (it1e.end ()); 00700 #else 00701 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00702 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00703 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00704 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00705 #endif 00706 difference_type it2_size (it2_end - it2); 00707 difference_type it2e_size (it2e_end - it2e); 00708 difference_type diff2 (0); 00709 if (it2_size > 0 && it2e_size > 0) { 00710 diff2 = it2.index2 () - it2e.index2 (); 00711 difference_type size2 = (std::min) (diff2, it2e_size); 00712 if (size2 > 0) { 00713 it2e += size2; 00714 it2e_size -= size2; 00715 diff2 -= size2; 00716 } 00717 size2 = (std::min) (- diff2, it2_size); 00718 if (size2 > 0) { 00719 it2_size -= size2; 00720 if (!functor_type::computed) { 00721 while (-- size2 >= 0) // zeroing 00722 functor_type::apply (*it2, expr_value_type/*zero*/()), ++ it2; 00723 } else { 00724 it2 += size2; 00725 } 00726 diff2 += size2; 00727 } 00728 } 00729 difference_type size2 ((std::min) (it2_size, it2e_size)); 00730 it2_size -= size2; 00731 it2e_size -= size2; 00732 while (-- size2 >= 0) 00733 functor_type::apply (*it2, *it2e), ++ it2, ++ it2e; 00734 size2 = it2_size; 00735 if (!functor_type::computed) { 00736 while (-- size2 >= 0) // zeroing 00737 functor_type::apply (*it2, expr_value_type/*zero*/()), ++ it2; 00738 } else { 00739 it2 += size2; 00740 } 00741 ++ it1, ++ it1e; 00742 } 00743 size1 = it1_size; 00744 if (!functor_type::computed) { 00745 while (-- size1 >= 0) { // zeroing 00746 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00747 typename M::iterator2 it2 (it1.begin ()); 00748 typename M::iterator2 it2_end (it1.end ()); 00749 #else 00750 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00751 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00752 #endif 00753 difference_type size2 (it2_end - it2); 00754 while (-- size2 >= 0) 00755 functor_type::apply (*it2, expr_value_type/*zero*/()), ++ it2; 00756 ++ it1; 00757 } 00758 } else { 00759 it1 += size1; 00760 } 00761 #if BOOST_UBLAS_TYPE_CHECK 00762 if (! disable_type_check<bool>::value) 00763 BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ()); 00764 #endif 00765 } 00766 // Packed (proxy) column major case 00767 template<template <class T1, class T2> class F, class R, class M, class E> 00768 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00769 void matrix_assign (M &m, const matrix_expression<E> &e, packed_proxy_tag, column_major_tag) { 00770 typedef typename matrix_traits<E>::value_type expr_value_type; 00771 typedef F<typename M::iterator1::reference, expr_value_type> functor_type; 00772 // R unnecessary, make_conformant not required 00773 typedef typename M::difference_type difference_type; 00774 typedef typename M::value_type value_type; 00775 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00776 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00777 #if BOOST_UBLAS_TYPE_CHECK 00778 matrix<value_type, column_major> cm (m.size1 (), m.size2 ()); 00779 indexing_matrix_assign<scalar_assign> (cm, m, column_major_tag ()); 00780 indexing_matrix_assign<F> (cm, e, column_major_tag ()); 00781 #endif 00782 typename M::iterator2 it2 (m.begin2 ()); 00783 typename M::iterator2 it2_end (m.end2 ()); 00784 typename E::const_iterator2 it2e (e ().begin2 ()); 00785 typename E::const_iterator2 it2e_end (e ().end2 ()); 00786 difference_type it2_size (it2_end - it2); 00787 difference_type it2e_size (it2e_end - it2e); 00788 difference_type diff2 (0); 00789 if (it2_size > 0 && it2e_size > 0) 00790 diff2 = it2.index2 () - it2e.index2 (); 00791 if (diff2 != 0) { 00792 difference_type size2 = (std::min) (diff2, it2e_size); 00793 if (size2 > 0) { 00794 it2e += size2; 00795 it2e_size -= size2; 00796 diff2 -= size2; 00797 } 00798 size2 = (std::min) (- diff2, it2_size); 00799 if (size2 > 0) { 00800 it2_size -= size2; 00801 if (!functor_type::computed) { 00802 while (-- size2 >= 0) { // zeroing 00803 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00804 typename M::iterator1 it1 (it2.begin ()); 00805 typename M::iterator1 it1_end (it2.end ()); 00806 #else 00807 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00808 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 00809 #endif 00810 difference_type size1 (it1_end - it1); 00811 while (-- size1 >= 0) 00812 functor_type::apply (*it1, expr_value_type/*zero*/()), ++ it1; 00813 ++ it2; 00814 } 00815 } else { 00816 it2 += size2; 00817 } 00818 diff2 += size2; 00819 } 00820 } 00821 difference_type size2 ((std::min) (it2_size, it2e_size)); 00822 it2_size -= size2; 00823 it2e_size -= size2; 00824 while (-- size2 >= 0) { 00825 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00826 typename M::iterator1 it1 (it2.begin ()); 00827 typename M::iterator1 it1_end (it2.end ()); 00828 typename E::const_iterator1 it1e (it2e.begin ()); 00829 typename E::const_iterator1 it1e_end (it2e.end ()); 00830 #else 00831 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00832 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 00833 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00834 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 00835 #endif 00836 difference_type it1_size (it1_end - it1); 00837 difference_type it1e_size (it1e_end - it1e); 00838 difference_type diff1 (0); 00839 if (it1_size > 0 && it1e_size > 0) { 00840 diff1 = it1.index1 () - it1e.index1 (); 00841 difference_type size1 = (std::min) (diff1, it1e_size); 00842 if (size1 > 0) { 00843 it1e += size1; 00844 it1e_size -= size1; 00845 diff1 -= size1; 00846 } 00847 size1 = (std::min) (- diff1, it1_size); 00848 if (size1 > 0) { 00849 it1_size -= size1; 00850 if (!functor_type::computed) { 00851 while (-- size1 >= 0) // zeroing 00852 functor_type::apply (*it1, expr_value_type/*zero*/()), ++ it1; 00853 } else { 00854 it1 += size1; 00855 } 00856 diff1 += size1; 00857 } 00858 } 00859 difference_type size1 ((std::min) (it1_size, it1e_size)); 00860 it1_size -= size1; 00861 it1e_size -= size1; 00862 while (-- size1 >= 0) 00863 functor_type::apply (*it1, *it1e), ++ it1, ++ it1e; 00864 size1 = it1_size; 00865 if (!functor_type::computed) { 00866 while (-- size1 >= 0) // zeroing 00867 functor_type::apply (*it1, expr_value_type/*zero*/()), ++ it1; 00868 } else { 00869 it1 += size1; 00870 } 00871 ++ it2, ++ it2e; 00872 } 00873 size2 = it2_size; 00874 if (!functor_type::computed) { 00875 while (-- size2 >= 0) { // zeroing 00876 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00877 typename M::iterator1 it1 (it2.begin ()); 00878 typename M::iterator1 it1_end (it2.end ()); 00879 #else 00880 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 00881 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 00882 #endif 00883 difference_type size1 (it1_end - it1); 00884 while (-- size1 >= 0) 00885 functor_type::apply (*it1, expr_value_type/*zero*/()), ++ it1; 00886 ++ it2; 00887 } 00888 } else { 00889 it2 += size2; 00890 } 00891 #if BOOST_UBLAS_TYPE_CHECK 00892 if (! disable_type_check<bool>::value) 00893 BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ()); 00894 #endif 00895 } 00896 // Sparse row major case 00897 template<template <class T1, class T2> class F, class R, class M, class E> 00898 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00899 void matrix_assign (M &m, const matrix_expression<E> &e, sparse_tag, row_major_tag) { 00900 typedef F<typename M::iterator2::reference, typename E::value_type> functor_type; 00901 // R unnecessary, make_conformant not required 00902 BOOST_STATIC_ASSERT ((!functor_type::computed)); 00903 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00904 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00905 typedef typename M::value_type value_type; 00906 // Sparse type has no numeric constraints to check 00907 00908 m.clear (); 00909 typename E::const_iterator1 it1e (e ().begin1 ()); 00910 typename E::const_iterator1 it1e_end (e ().end1 ()); 00911 while (it1e != it1e_end) { 00912 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00913 typename E::const_iterator2 it2e (it1e.begin ()); 00914 typename E::const_iterator2 it2e_end (it1e.end ()); 00915 #else 00916 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00917 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00918 #endif 00919 while (it2e != it2e_end) { 00920 value_type t (*it2e); 00921 if (t != value_type/*zero*/()) 00922 m.insert_element (it2e.index1 (), it2e.index2 (), t); 00923 ++ it2e; 00924 } 00925 ++ it1e; 00926 } 00927 } 00928 // Sparse column major case 00929 template<template <class T1, class T2> class F, class R, class M, class E> 00930 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00931 void matrix_assign (M &m, const matrix_expression<E> &e, sparse_tag, column_major_tag) { 00932 typedef F<typename M::iterator1::reference, typename E::value_type> functor_type; 00933 // R unnecessary, make_conformant not required 00934 BOOST_STATIC_ASSERT ((!functor_type::computed)); 00935 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00936 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00937 typedef typename M::value_type value_type; 00938 // Sparse type has no numeric constraints to check 00939 00940 m.clear (); 00941 typename E::const_iterator2 it2e (e ().begin2 ()); 00942 typename E::const_iterator2 it2e_end (e ().end2 ()); 00943 while (it2e != it2e_end) { 00944 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00945 typename E::const_iterator1 it1e (it2e.begin ()); 00946 typename E::const_iterator1 it1e_end (it2e.end ()); 00947 #else 00948 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 00949 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 00950 #endif 00951 while (it1e != it1e_end) { 00952 value_type t (*it1e); 00953 if (t != value_type/*zero*/()) 00954 m.insert_element (it1e.index1 (), it1e.index2 (), t); 00955 ++ it1e; 00956 } 00957 ++ it2e; 00958 } 00959 } 00960 // Sparse proxy or functional row major case 00961 template<template <class T1, class T2> class F, class R, class M, class E> 00962 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 00963 void matrix_assign (M &m, const matrix_expression<E> &e, sparse_proxy_tag, row_major_tag) { 00964 typedef typename matrix_traits<E>::value_type expr_value_type; 00965 typedef F<typename M::iterator2::reference, expr_value_type> functor_type; 00966 typedef R conformant_restrict_type; 00967 typedef typename M::size_type size_type; 00968 typedef typename M::difference_type difference_type; 00969 typedef typename M::value_type value_type; 00970 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 00971 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 00972 #if BOOST_UBLAS_TYPE_CHECK 00973 matrix<value_type, row_major> cm (m.size1 (), m.size2 ()); 00974 indexing_matrix_assign<scalar_assign> (cm, m, row_major_tag ()); 00975 indexing_matrix_assign<F> (cm, e, row_major_tag ()); 00976 #endif 00977 detail::make_conformant (m, e, row_major_tag (), conformant_restrict_type ()); 00978 00979 typename M::iterator1 it1 (m.begin1 ()); 00980 typename M::iterator1 it1_end (m.end1 ()); 00981 typename E::const_iterator1 it1e (e ().begin1 ()); 00982 typename E::const_iterator1 it1e_end (e ().end1 ()); 00983 while (it1 != it1_end && it1e != it1e_end) { 00984 difference_type compare = it1.index1 () - it1e.index1 (); 00985 if (compare == 0) { 00986 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 00987 typename M::iterator2 it2 (it1.begin ()); 00988 typename M::iterator2 it2_end (it1.end ()); 00989 typename E::const_iterator2 it2e (it1e.begin ()); 00990 typename E::const_iterator2 it2e_end (it1e.end ()); 00991 #else 00992 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 00993 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 00994 typename E::const_iterator2 it2e (begin (it1e, iterator1_tag ())); 00995 typename E::const_iterator2 it2e_end (end (it1e, iterator1_tag ())); 00996 #endif 00997 if (it2 != it2_end && it2e != it2e_end) { 00998 size_type it2_index = it2.index2 (), it2e_index = it2e.index2 (); 00999 while (true) { 01000 difference_type compare = it2_index - it2e_index; 01001 if (compare == 0) { 01002 functor_type::apply (*it2, *it2e); 01003 ++ it2, ++ it2e; 01004 if (it2 != it2_end && it2e != it2e_end) { 01005 it2_index = it2.index2 (); 01006 it2e_index = it2e.index2 (); 01007 } else 01008 break; 01009 } else if (compare < 0) { 01010 if (!functor_type::computed) { 01011 functor_type::apply (*it2, expr_value_type/*zero*/()); 01012 ++ it2; 01013 } else 01014 increment (it2, it2_end, - compare); 01015 if (it2 != it2_end) 01016 it2_index = it2.index2 (); 01017 else 01018 break; 01019 } else if (compare > 0) { 01020 increment (it2e, it2e_end, compare); 01021 if (it2e != it2e_end) 01022 it2e_index = it2e.index2 (); 01023 else 01024 break; 01025 } 01026 } 01027 } 01028 if (!functor_type::computed) { 01029 while (it2 != it2_end) { // zeroing 01030 functor_type::apply (*it2, expr_value_type/*zero*/()); 01031 ++ it2; 01032 } 01033 } else { 01034 it2 = it2_end; 01035 } 01036 ++ it1, ++ it1e; 01037 } else if (compare < 0) { 01038 if (!functor_type::computed) { 01039 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01040 typename M::iterator2 it2 (it1.begin ()); 01041 typename M::iterator2 it2_end (it1.end ()); 01042 #else 01043 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01044 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 01045 #endif 01046 while (it2 != it2_end) { // zeroing 01047 functor_type::apply (*it2, expr_value_type/*zero*/()); 01048 ++ it2; 01049 } 01050 ++ it1; 01051 } else { 01052 increment (it1, it1_end, - compare); 01053 } 01054 } else if (compare > 0) { 01055 increment (it1e, it1e_end, compare); 01056 } 01057 } 01058 if (!functor_type::computed) { 01059 while (it1 != it1_end) { 01060 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01061 typename M::iterator2 it2 (it1.begin ()); 01062 typename M::iterator2 it2_end (it1.end ()); 01063 #else 01064 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01065 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 01066 #endif 01067 while (it2 != it2_end) { // zeroing 01068 functor_type::apply (*it2, expr_value_type/*zero*/()); 01069 ++ it2; 01070 } 01071 ++ it1; 01072 } 01073 } else { 01074 it1 = it1_end; 01075 } 01076 #if BOOST_UBLAS_TYPE_CHECK 01077 if (! disable_type_check<bool>::value) 01078 BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ()); 01079 #endif 01080 } 01081 // Sparse proxy or functional column major case 01082 template<template <class T1, class T2> class F, class R, class M, class E> 01083 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01084 void matrix_assign (M &m, const matrix_expression<E> &e, sparse_proxy_tag, column_major_tag) { 01085 typedef typename matrix_traits<E>::value_type expr_value_type; 01086 typedef F<typename M::iterator1::reference, expr_value_type> functor_type; 01087 typedef R conformant_restrict_type; 01088 typedef typename M::size_type size_type; 01089 typedef typename M::difference_type difference_type; 01090 typedef typename M::value_type value_type; 01091 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 01092 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 01093 #if BOOST_UBLAS_TYPE_CHECK 01094 matrix<value_type, column_major> cm (m.size1 (), m.size2 ()); 01095 indexing_matrix_assign<scalar_assign> (cm, m, column_major_tag ()); 01096 indexing_matrix_assign<F> (cm, e, column_major_tag ()); 01097 #endif 01098 detail::make_conformant (m, e, column_major_tag (), conformant_restrict_type ()); 01099 01100 typename M::iterator2 it2 (m.begin2 ()); 01101 typename M::iterator2 it2_end (m.end2 ()); 01102 typename E::const_iterator2 it2e (e ().begin2 ()); 01103 typename E::const_iterator2 it2e_end (e ().end2 ()); 01104 while (it2 != it2_end && it2e != it2e_end) { 01105 difference_type compare = it2.index2 () - it2e.index2 (); 01106 if (compare == 0) { 01107 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01108 typename M::iterator1 it1 (it2.begin ()); 01109 typename M::iterator1 it1_end (it2.end ()); 01110 typename E::const_iterator1 it1e (it2e.begin ()); 01111 typename E::const_iterator1 it1e_end (it2e.end ()); 01112 #else 01113 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01114 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01115 typename E::const_iterator1 it1e (begin (it2e, iterator2_tag ())); 01116 typename E::const_iterator1 it1e_end (end (it2e, iterator2_tag ())); 01117 #endif 01118 if (it1 != it1_end && it1e != it1e_end) { 01119 size_type it1_index = it1.index1 (), it1e_index = it1e.index1 (); 01120 while (true) { 01121 difference_type compare = it1_index - it1e_index; 01122 if (compare == 0) { 01123 functor_type::apply (*it1, *it1e); 01124 ++ it1, ++ it1e; 01125 if (it1 != it1_end && it1e != it1e_end) { 01126 it1_index = it1.index1 (); 01127 it1e_index = it1e.index1 (); 01128 } else 01129 break; 01130 } else if (compare < 0) { 01131 if (!functor_type::computed) { 01132 functor_type::apply (*it1, expr_value_type/*zero*/()); // zeroing 01133 ++ it1; 01134 } else 01135 increment (it1, it1_end, - compare); 01136 if (it1 != it1_end) 01137 it1_index = it1.index1 (); 01138 else 01139 break; 01140 } else if (compare > 0) { 01141 increment (it1e, it1e_end, compare); 01142 if (it1e != it1e_end) 01143 it1e_index = it1e.index1 (); 01144 else 01145 break; 01146 } 01147 } 01148 } 01149 if (!functor_type::computed) { 01150 while (it1 != it1_end) { // zeroing 01151 functor_type::apply (*it1, expr_value_type/*zero*/()); 01152 ++ it1; 01153 } 01154 } else { 01155 it1 = it1_end; 01156 } 01157 ++ it2, ++ it2e; 01158 } else if (compare < 0) { 01159 if (!functor_type::computed) { 01160 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01161 typename M::iterator1 it1 (it2.begin ()); 01162 typename M::iterator1 it1_end (it2.end ()); 01163 #else 01164 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01165 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01166 #endif 01167 while (it1 != it1_end) { // zeroing 01168 functor_type::apply (*it1, expr_value_type/*zero*/()); 01169 ++ it1; 01170 } 01171 ++ it2; 01172 } else { 01173 increment (it2, it2_end, - compare); 01174 } 01175 } else if (compare > 0) { 01176 increment (it2e, it2e_end, compare); 01177 } 01178 } 01179 if (!functor_type::computed) { 01180 while (it2 != it2_end) { 01181 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01182 typename M::iterator1 it1 (it2.begin ()); 01183 typename M::iterator1 it1_end (it2.end ()); 01184 #else 01185 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01186 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01187 #endif 01188 while (it1 != it1_end) { // zeroing 01189 functor_type::apply (*it1, expr_value_type/*zero*/()); 01190 ++ it1; 01191 } 01192 ++ it2; 01193 } 01194 } else { 01195 it2 = it2_end; 01196 } 01197 #if BOOST_UBLAS_TYPE_CHECK 01198 if (! disable_type_check<bool>::value) 01199 BOOST_UBLAS_CHECK (detail::expression_type_check (m, cm), external_logic ()); 01200 #endif 01201 } 01202 01203 // Dispatcher 01204 template<template <class T1, class T2> class F, class M, class E> 01205 BOOST_UBLAS_INLINE 01206 void matrix_assign (M &m, const matrix_expression<E> &e) { 01207 typedef typename matrix_assign_traits<typename M::storage_category, 01208 F<typename M::reference, typename E::value_type>::computed, 01209 typename E::const_iterator1::iterator_category, 01210 typename E::const_iterator2::iterator_category>::storage_category storage_category; 01211 // give preference to matrix M's orientation if known 01212 typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>, 01213 typename E::orientation_category , 01214 typename M::orientation_category >::type orientation_category; 01215 typedef basic_full<typename M::size_type> unrestricted; 01216 matrix_assign<F, unrestricted> (m, e, storage_category (), orientation_category ()); 01217 } 01218 template<template <class T1, class T2> class F, class R, class M, class E> 01219 BOOST_UBLAS_INLINE 01220 void matrix_assign (M &m, const matrix_expression<E> &e) { 01221 typedef R conformant_restrict_type; 01222 typedef typename matrix_assign_traits<typename M::storage_category, 01223 F<typename M::reference, typename E::value_type>::computed, 01224 typename E::const_iterator1::iterator_category, 01225 typename E::const_iterator2::iterator_category>::storage_category storage_category; 01226 // give preference to matrix M's orientation if known 01227 typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>, 01228 typename E::orientation_category , 01229 typename M::orientation_category >::type orientation_category; 01230 matrix_assign<F, conformant_restrict_type> (m, e, storage_category (), orientation_category ()); 01231 } 01232 01233 template<class SC, class RI1, class RI2> 01234 struct matrix_swap_traits { 01235 typedef SC storage_category; 01236 }; 01237 01238 template<> 01239 struct matrix_swap_traits<dense_proxy_tag, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 01240 typedef sparse_proxy_tag storage_category; 01241 }; 01242 01243 template<> 01244 struct matrix_swap_traits<packed_proxy_tag, sparse_bidirectional_iterator_tag, sparse_bidirectional_iterator_tag> { 01245 typedef sparse_proxy_tag storage_category; 01246 }; 01247 01248 // Dense (proxy) row major case 01249 template<template <class T1, class T2> class F, class R, class M, class E> 01250 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01251 void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, row_major_tag) { 01252 typedef F<typename M::iterator2::reference, typename E::reference> functor_type; 01253 // R unnecessary, make_conformant not required 01254 typedef typename M::size_type size_type; 01255 typedef typename M::difference_type difference_type; 01256 typename M::iterator1 it1 (m.begin1 ()); 01257 typename E::iterator1 it1e (e ().begin1 ()); 01258 difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (e ().end1 () - it1e))); 01259 while (-- size1 >= 0) { 01260 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01261 typename M::iterator2 it2 (it1.begin ()); 01262 typename E::iterator2 it2e (it1e.begin ()); 01263 difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (it1e.end () - it2e))); 01264 #else 01265 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01266 typename E::iterator2 it2e (begin (it1e, iterator1_tag ())); 01267 difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (end (it1e, iterator1_tag ()) - it2e))); 01268 #endif 01269 while (-- size2 >= 0) 01270 functor_type::apply (*it2, *it2e), ++ it2, ++ it2e; 01271 ++ it1, ++ it1e; 01272 } 01273 } 01274 // Dense (proxy) column major case 01275 template<template <class T1, class T2> class F, class R, class M, class E> 01276 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01277 void matrix_swap (M &m, matrix_expression<E> &e, dense_proxy_tag, column_major_tag) { 01278 typedef F<typename M::iterator1::reference, typename E::reference> functor_type; 01279 // R unnecessary, make_conformant not required 01280 typedef typename M::size_type size_type; 01281 typedef typename M::difference_type difference_type; 01282 typename M::iterator2 it2 (m.begin2 ()); 01283 typename E::iterator2 it2e (e ().begin2 ()); 01284 difference_type size2 (BOOST_UBLAS_SAME (m.size2 (), size_type (e ().end2 () - it2e))); 01285 while (-- size2 >= 0) { 01286 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01287 typename M::iterator1 it1 (it2.begin ()); 01288 typename E::iterator1 it1e (it2e.begin ()); 01289 difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (it2e.end () - it1e))); 01290 #else 01291 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01292 typename E::iterator1 it1e (begin (it2e, iterator2_tag ())); 01293 difference_type size1 (BOOST_UBLAS_SAME (m.size1 (), size_type (end (it2e, iterator2_tag ()) - it1e))); 01294 #endif 01295 while (-- size1 >= 0) 01296 functor_type::apply (*it1, *it1e), ++ it1, ++ it1e; 01297 ++ it2, ++ it2e; 01298 } 01299 } 01300 // Packed (proxy) row major case 01301 template<template <class T1, class T2> class F, class R, class M, class E> 01302 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01303 void matrix_swap (M &m, matrix_expression<E> &e, packed_proxy_tag, row_major_tag) { 01304 typedef F<typename M::iterator2::reference, typename E::reference> functor_type; 01305 // R unnecessary, make_conformant not required 01306 typedef typename M::size_type size_type; 01307 typedef typename M::difference_type difference_type; 01308 typename M::iterator1 it1 (m.begin1 ()); 01309 typename E::iterator1 it1e (e ().begin1 ()); 01310 difference_type size1 (BOOST_UBLAS_SAME (m.end1 () - it1, e ().end1 () - it1e)); 01311 while (-- size1 >= 0) { 01312 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01313 typename M::iterator2 it2 (it1.begin ()); 01314 typename E::iterator2 it2e (it1e.begin ()); 01315 difference_type size2 (BOOST_UBLAS_SAME (it1.end () - it2, it1e.end () - it2e)); 01316 #else 01317 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01318 typename E::iterator2 it2e (begin (it1e, iterator1_tag ())); 01319 difference_type size2 (BOOST_UBLAS_SAME (end (it1, iterator1_tag ()) - it2, end (it1e, iterator1_tag ()) - it2e)); 01320 #endif 01321 while (-- size2 >= 0) 01322 functor_type::apply (*it2, *it2e), ++ it2, ++ it2e; 01323 ++ it1, ++ it1e; 01324 } 01325 } 01326 // Packed (proxy) column major case 01327 template<template <class T1, class T2> class F, class R, class M, class E> 01328 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01329 void matrix_swap (M &m, matrix_expression<E> &e, packed_proxy_tag, column_major_tag) { 01330 typedef F<typename M::iterator1::reference, typename E::reference> functor_type; 01331 // R unnecessary, make_conformant not required 01332 typedef typename M::size_type size_type; 01333 typedef typename M::difference_type difference_type; 01334 typename M::iterator2 it2 (m.begin2 ()); 01335 typename E::iterator2 it2e (e ().begin2 ()); 01336 difference_type size2 (BOOST_UBLAS_SAME (m.end2 () - it2, e ().end2 () - it2e)); 01337 while (-- size2 >= 0) { 01338 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01339 typename M::iterator1 it1 (it2.begin ()); 01340 typename E::iterator1 it1e (it2e.begin ()); 01341 difference_type size1 (BOOST_UBLAS_SAME (it2.end () - it1, it2e.end () - it1e)); 01342 #else 01343 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01344 typename E::iterator1 it1e (begin (it2e, iterator2_tag ())); 01345 difference_type size1 (BOOST_UBLAS_SAME (end (it2, iterator2_tag ()) - it1, end (it2e, iterator2_tag ()) - it1e)); 01346 #endif 01347 while (-- size1 >= 0) 01348 functor_type::apply (*it1, *it1e), ++ it1, ++ it1e; 01349 ++ it2, ++ it2e; 01350 } 01351 } 01352 // Sparse (proxy) row major case 01353 template<template <class T1, class T2> class F, class R, class M, class E> 01354 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01355 void matrix_swap (M &m, matrix_expression<E> &e, sparse_proxy_tag, row_major_tag) { 01356 typedef F<typename M::iterator2::reference, typename E::reference> functor_type; 01357 typedef R conformant_restrict_type; 01358 typedef typename M::size_type size_type; 01359 typedef typename M::difference_type difference_type; 01360 typedef typename M::value_type value_type; 01361 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 01362 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 01363 01364 detail::make_conformant (m, e, row_major_tag (), conformant_restrict_type ()); 01365 // FIXME should be a seperate restriction for E 01366 detail::make_conformant (e (), m, row_major_tag (), conformant_restrict_type ()); 01367 01368 typename M::iterator1 it1 (m.begin1 ()); 01369 typename M::iterator1 it1_end (m.end1 ()); 01370 typename E::iterator1 it1e (e ().begin1 ()); 01371 typename E::iterator1 it1e_end (e ().end1 ()); 01372 while (it1 != it1_end && it1e != it1e_end) { 01373 difference_type compare = it1.index1 () - it1e.index1 (); 01374 if (compare == 0) { 01375 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01376 typename M::iterator2 it2 (it1.begin ()); 01377 typename M::iterator2 it2_end (it1.end ()); 01378 typename E::iterator2 it2e (it1e.begin ()); 01379 typename E::iterator2 it2e_end (it1e.end ()); 01380 #else 01381 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01382 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 01383 typename E::iterator2 it2e (begin (it1e, iterator1_tag ())); 01384 typename E::iterator2 it2e_end (end (it1e, iterator1_tag ())); 01385 #endif 01386 if (it2 != it2_end && it2e != it2e_end) { 01387 size_type it2_index = it2.index2 (), it2e_index = it2e.index2 (); 01388 while (true) { 01389 difference_type compare = it2_index - it2e_index; 01390 if (compare == 0) { 01391 functor_type::apply (*it2, *it2e); 01392 ++ it2, ++ it2e; 01393 if (it2 != it2_end && it2e != it2e_end) { 01394 it2_index = it2.index2 (); 01395 it2e_index = it2e.index2 (); 01396 } else 01397 break; 01398 } else if (compare < 0) { 01399 increment (it2, it2_end, - compare); 01400 if (it2 != it2_end) 01401 it2_index = it2.index2 (); 01402 else 01403 break; 01404 } else if (compare > 0) { 01405 increment (it2e, it2e_end, compare); 01406 if (it2e != it2e_end) 01407 it2e_index = it2e.index2 (); 01408 else 01409 break; 01410 } 01411 } 01412 } 01413 #if BOOST_UBLAS_TYPE_CHECK 01414 increment (it2e, it2e_end); 01415 increment (it2, it2_end); 01416 #endif 01417 ++ it1, ++ it1e; 01418 } else if (compare < 0) { 01419 #if BOOST_UBLAS_TYPE_CHECK 01420 while (it1.index1 () < it1e.index1 ()) { 01421 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01422 typename M::iterator2 it2 (it1.begin ()); 01423 typename M::iterator2 it2_end (it1.end ()); 01424 #else 01425 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01426 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 01427 #endif 01428 increment (it2, it2_end); 01429 ++ it1; 01430 } 01431 #else 01432 increment (it1, it1_end, - compare); 01433 #endif 01434 } else if (compare > 0) { 01435 #if BOOST_UBLAS_TYPE_CHECK 01436 while (it1e.index1 () < it1.index1 ()) { 01437 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01438 typename E::iterator2 it2e (it1e.begin ()); 01439 typename E::iterator2 it2e_end (it1e.end ()); 01440 #else 01441 typename E::iterator2 it2e (begin (it1e, iterator1_tag ())); 01442 typename E::iterator2 it2e_end (end (it1e, iterator1_tag ())); 01443 #endif 01444 increment (it2e, it2e_end); 01445 ++ it1e; 01446 } 01447 #else 01448 increment (it1e, it1e_end, compare); 01449 #endif 01450 } 01451 } 01452 #if BOOST_UBLAS_TYPE_CHECK 01453 while (it1e != it1e_end) { 01454 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01455 typename E::iterator2 it2e (it1e.begin ()); 01456 typename E::iterator2 it2e_end (it1e.end ()); 01457 #else 01458 typename E::iterator2 it2e (begin (it1e, iterator1_tag ())); 01459 typename E::iterator2 it2e_end (end (it1e, iterator1_tag ())); 01460 #endif 01461 increment (it2e, it2e_end); 01462 ++ it1e; 01463 } 01464 while (it1 != it1_end) { 01465 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01466 typename M::iterator2 it2 (it1.begin ()); 01467 typename M::iterator2 it2_end (it1.end ()); 01468 #else 01469 typename M::iterator2 it2 (begin (it1, iterator1_tag ())); 01470 typename M::iterator2 it2_end (end (it1, iterator1_tag ())); 01471 #endif 01472 increment (it2, it2_end); 01473 ++ it1; 01474 } 01475 #endif 01476 } 01477 // Sparse (proxy) column major case 01478 template<template <class T1, class T2> class F, class R, class M, class E> 01479 // BOOST_UBLAS_INLINE This function seems to be big. So we do not let the compiler inline it. 01480 void matrix_swap (M &m, matrix_expression<E> &e, sparse_proxy_tag, column_major_tag) { 01481 typedef F<typename M::iterator1::reference, typename E::reference> functor_type; 01482 typedef R conformant_restrict_type; 01483 typedef typename M::size_type size_type; 01484 typedef typename M::difference_type difference_type; 01485 typedef typename M::value_type value_type; 01486 BOOST_UBLAS_CHECK (m.size1 () == e ().size1 (), bad_size ()); 01487 BOOST_UBLAS_CHECK (m.size2 () == e ().size2 (), bad_size ()); 01488 01489 detail::make_conformant (m, e, column_major_tag (), conformant_restrict_type ()); 01490 // FIXME should be a seperate restriction for E 01491 detail::make_conformant (e (), m, column_major_tag (), conformant_restrict_type ()); 01492 01493 typename M::iterator2 it2 (m.begin2 ()); 01494 typename M::iterator2 it2_end (m.end2 ()); 01495 typename E::iterator2 it2e (e ().begin2 ()); 01496 typename E::iterator2 it2e_end (e ().end2 ()); 01497 while (it2 != it2_end && it2e != it2e_end) { 01498 difference_type compare = it2.index2 () - it2e.index2 (); 01499 if (compare == 0) { 01500 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01501 typename M::iterator1 it1 (it2.begin ()); 01502 typename M::iterator1 it1_end (it2.end ()); 01503 typename E::iterator1 it1e (it2e.begin ()); 01504 typename E::iterator1 it1e_end (it2e.end ()); 01505 #else 01506 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01507 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01508 typename E::iterator1 it1e (begin (it2e, iterator2_tag ())); 01509 typename E::iterator1 it1e_end (end (it2e, iterator2_tag ())); 01510 #endif 01511 if (it1 != it1_end && it1e != it1e_end) { 01512 size_type it1_index = it1.index1 (), it1e_index = it1e.index1 (); 01513 while (true) { 01514 difference_type compare = it1_index - it1e_index; 01515 if (compare == 0) { 01516 functor_type::apply (*it1, *it1e); 01517 ++ it1, ++ it1e; 01518 if (it1 != it1_end && it1e != it1e_end) { 01519 it1_index = it1.index1 (); 01520 it1e_index = it1e.index1 (); 01521 } else 01522 break; 01523 } else if (compare < 0) { 01524 increment (it1, it1_end, - compare); 01525 if (it1 != it1_end) 01526 it1_index = it1.index1 (); 01527 else 01528 break; 01529 } else if (compare > 0) { 01530 increment (it1e, it1e_end, compare); 01531 if (it1e != it1e_end) 01532 it1e_index = it1e.index1 (); 01533 else 01534 break; 01535 } 01536 } 01537 } 01538 #if BOOST_UBLAS_TYPE_CHECK 01539 increment (it1e, it1e_end); 01540 increment (it1, it1_end); 01541 #endif 01542 ++ it2, ++ it2e; 01543 } else if (compare < 0) { 01544 #if BOOST_UBLAS_TYPE_CHECK 01545 while (it2.index2 () < it2e.index2 ()) { 01546 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01547 typename M::iterator1 it1 (it2.begin ()); 01548 typename M::iterator1 it1_end (it2.end ()); 01549 #else 01550 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01551 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01552 #endif 01553 increment (it1, it1_end); 01554 ++ it2; 01555 } 01556 #else 01557 increment (it2, it2_end, - compare); 01558 #endif 01559 } else if (compare > 0) { 01560 #if BOOST_UBLAS_TYPE_CHECK 01561 while (it2e.index2 () < it2.index2 ()) { 01562 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01563 typename E::iterator1 it1e (it2e.begin ()); 01564 typename E::iterator1 it1e_end (it2e.end ()); 01565 #else 01566 typename E::iterator1 it1e (begin (it2e, iterator2_tag ())); 01567 typename E::iterator1 it1e_end (end (it2e, iterator2_tag ())); 01568 #endif 01569 increment (it1e, it1e_end); 01570 ++ it2e; 01571 } 01572 #else 01573 increment (it2e, it2e_end, compare); 01574 #endif 01575 } 01576 } 01577 #if BOOST_UBLAS_TYPE_CHECK 01578 while (it2e != it2e_end) { 01579 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01580 typename E::iterator1 it1e (it2e.begin ()); 01581 typename E::iterator1 it1e_end (it2e.end ()); 01582 #else 01583 typename E::iterator1 it1e (begin (it2e, iterator2_tag ())); 01584 typename E::iterator1 it1e_end (end (it2e, iterator2_tag ())); 01585 #endif 01586 increment (it1e, it1e_end); 01587 ++ it2e; 01588 } 01589 while (it2 != it2_end) { 01590 #ifndef BOOST_UBLAS_NO_NESTED_CLASS_RELATION 01591 typename M::iterator1 it1 (it2.begin ()); 01592 typename M::iterator1 it1_end (it2.end ()); 01593 #else 01594 typename M::iterator1 it1 (begin (it2, iterator2_tag ())); 01595 typename M::iterator1 it1_end (end (it2, iterator2_tag ())); 01596 #endif 01597 increment (it1, it1_end); 01598 ++ it2; 01599 } 01600 #endif 01601 } 01602 01603 // Dispatcher 01604 template<template <class T1, class T2> class F, class M, class E> 01605 BOOST_UBLAS_INLINE 01606 void matrix_swap (M &m, matrix_expression<E> &e) { 01607 typedef typename matrix_swap_traits<typename M::storage_category, 01608 typename E::const_iterator1::iterator_category, 01609 typename E::const_iterator2::iterator_category>::storage_category storage_category; 01610 // give preference to matrix M's orientation if known 01611 typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>, 01612 typename E::orientation_category , 01613 typename M::orientation_category >::type orientation_category; 01614 typedef basic_full<typename M::size_type> unrestricted; 01615 matrix_swap<F, unrestricted> (m, e, storage_category (), orientation_category ()); 01616 } 01617 template<template <class T1, class T2> class F, class R, class M, class E> 01618 BOOST_UBLAS_INLINE 01619 void matrix_swap (M &m, matrix_expression<E> &e) { 01620 typedef R conformant_restrict_type; 01621 typedef typename matrix_swap_traits<typename M::storage_category, 01622 typename E::const_iterator1::iterator_category, 01623 typename E::const_iterator2::iterator_category>::storage_category storage_category; 01624 // give preference to matrix M's orientation if known 01625 typedef typename boost::mpl::if_<boost::is_same<typename M::orientation_category, unknown_orientation_tag>, 01626 typename E::orientation_category , 01627 typename M::orientation_category >::type orientation_category; 01628 matrix_swap<F, conformant_restrict_type> (m, e, storage_category (), orientation_category ()); 01629 } 01630 01631 }}} 01632 01633 #endif