2 // Copyright (c) 2000-2002
3 // Joerg Walter, Mathias Koch
5 // Distributed under the Boost Software License, Version 1.0. (See
6 // accompanying file LICENSE_1_0.txt or copy at
7 // http://www.boost.org/LICENSE_1_0.txt)
9 // The authors gratefully acknowledge the support of
10 // GeNeSys mbH & Co. KG in producing this work.
13 #ifndef _BOOST_UBLAS_VECTOR_EXPRESSION_
14 #define _BOOST_UBLAS_VECTOR_EXPRESSION_
16 #include <boost/numeric/ublas/expression_types.hpp>
19 // Expression templates based on ideas of Todd Veldhuizen and Geoffrey Furnish
20 // Iterators based on ideas of Jeremy Siek
22 // Classes that model the Vector Expression concept
24 namespace boost { namespace numeric { namespace ublas {
27 class vector_reference:
28 public vector_expression<vector_reference<E> > {
30 typedef vector_reference<E> self_type;
32 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
33 using vector_expression<vector_reference<E> >::operator ();
35 typedef typename E::size_type size_type;
36 typedef typename E::difference_type difference_type;
37 typedef typename E::value_type value_type;
38 typedef typename E::const_reference const_reference;
39 typedef typename boost::mpl::if_<boost::is_const<E>,
40 typename E::const_reference,
41 typename E::reference>::type reference;
42 typedef E referred_type;
43 typedef const self_type const_closure_type;
44 typedef self_type closure_type;
45 typedef typename E::storage_category storage_category;
47 // Construction and destruction
49 explicit vector_reference (referred_type &e):
54 size_type size () const {
55 return expression ().size ();
59 // Expression accessors - const correct
61 const referred_type &expression () const {
65 referred_type &expression () {
71 #ifndef BOOST_UBLAS_REFERENCE_CONST_MEMBER
73 const_reference operator () (size_type i) const {
74 return expression () (i);
77 reference operator () (size_type i) {
78 return expression () (i);
82 const_reference operator [] (size_type i) const {
83 return expression () [i];
86 reference operator [] (size_type i) {
87 return expression () [i];
91 reference operator () (size_type i) const {
92 return expression () (i);
96 reference operator [] (size_type i) const {
97 return expression () [i];
103 vector_reference &operator = (const vector_reference &v) {
104 expression ().operator = (v);
109 vector_reference &operator = (const vector_expression<AE> &ae) {
110 expression ().operator = (ae);
115 vector_reference &assign (const vector_expression<AE> &ae) {
116 expression ().assign (ae);
121 vector_reference &operator += (const vector_expression<AE> &ae) {
122 expression ().operator += (ae);
127 vector_reference &plus_assign (const vector_expression<AE> &ae) {
128 expression ().plus_assign (ae);
133 vector_reference &operator -= (const vector_expression<AE> &ae) {
134 expression ().operator -= (ae);
139 vector_reference &minus_assign (const vector_expression<AE> &ae) {
140 expression ().minus_assign (ae);
145 vector_reference &operator *= (const AT &at) {
146 expression ().operator *= (at);
151 vector_reference &operator /= (const AT &at) {
152 expression ().operator /= (at);
158 void swap (vector_reference &v) {
159 expression ().swap (v.expression ());
162 // Closure comparison
164 bool same_closure (const vector_reference &vr) const {
165 return &(*this).e_ == &vr.e_;
169 typedef typename E::const_iterator const_iterator;
170 typedef typename boost::mpl::if_<boost::is_const<E>,
171 typename E::const_iterator,
172 typename E::iterator>::type iterator;
176 const_iterator find (size_type i) const {
177 return expression ().find (i);
180 iterator find (size_type i) {
181 return expression ().find (i);
184 // Iterator is the iterator of the referenced expression.
187 const_iterator begin () const {
188 return expression ().begin ();
191 const_iterator cbegin () const {
195 const_iterator end () const {
196 return expression ().end ();
199 const_iterator cend () const {
205 return expression ().begin ();
209 return expression ().end ();
213 typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
214 typedef reverse_iterator_base<iterator> reverse_iterator;
217 const_reverse_iterator rbegin () const {
218 return const_reverse_iterator (end ());
221 const_reverse_iterator crbegin () const {
225 const_reverse_iterator rend () const {
226 return const_reverse_iterator (begin ());
229 const_reverse_iterator crend () const {
234 reverse_iterator rbegin () {
235 return reverse_iterator (end ());
238 reverse_iterator rend () {
239 return reverse_iterator (begin ());
247 template<class E, class F>
249 public vector_expression<vector_unary<E, F> > {
251 typedef F functor_type;
252 typedef typename boost::mpl::if_<boost::is_same<F, scalar_identity<typename E::value_type> >,
254 const E>::type expression_type;
255 typedef typename boost::mpl::if_<boost::is_const<expression_type>,
256 typename E::const_closure_type,
257 typename E::closure_type>::type expression_closure_type;
258 typedef vector_unary<E, F> self_type;
260 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
261 using vector_expression<vector_unary<E, F> >::operator ();
263 typedef typename E::size_type size_type;
264 typedef typename E::difference_type difference_type;
265 typedef typename F::result_type value_type;
266 typedef value_type const_reference;
267 typedef typename boost::mpl::if_<boost::is_same<F, scalar_identity<value_type> >,
268 typename E::reference,
269 value_type>::type reference;
270 typedef const self_type const_closure_type;
271 typedef self_type closure_type;
272 typedef unknown_storage_tag storage_category;
274 // Construction and destruction
276 // May be used as mutable expression.
277 explicit vector_unary (expression_type &e):
282 size_type size () const {
287 // Expression accessors
289 const expression_closure_type &expression () const {
296 const_reference operator () (size_type i) const {
297 return functor_type::apply (e_ (i));
300 reference operator () (size_type i) {
301 BOOST_STATIC_ASSERT ((boost::is_same<functor_type, scalar_identity<value_type > >::value));
306 const_reference operator [] (size_type i) const {
307 return functor_type::apply (e_ [i]);
310 reference operator [] (size_type i) {
311 BOOST_STATIC_ASSERT ((boost::is_same<functor_type, scalar_identity<value_type > >::value));
315 // Closure comparison
317 bool same_closure (const vector_unary &vu) const {
318 return (*this).expression ().same_closure (vu.expression ());
323 typedef typename E::const_iterator const_subiterator_type;
324 typedef const value_type *const_pointer;
327 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
328 typedef indexed_const_iterator<const_closure_type, typename const_subiterator_type::iterator_category> const_iterator;
329 typedef const_iterator iterator;
331 class const_iterator;
332 typedef const_iterator iterator;
337 const_iterator find (size_type i) const {
338 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
339 const_subiterator_type it (e_.find (i));
340 return const_iterator (*this, it.index ());
342 return const_iterator (*this, e_.find (i));
346 // Iterator enhances the iterator of the referenced expression
347 // with the unary functor.
349 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
350 class const_iterator:
351 public container_const_reference<vector_unary>,
352 public iterator_base_traits<typename E::const_iterator::iterator_category>::template
353 iterator_base<const_iterator, value_type>::type {
355 typedef typename E::const_iterator::iterator_category iterator_category;
356 typedef typename vector_unary::difference_type difference_type;
357 typedef typename vector_unary::value_type value_type;
358 typedef typename vector_unary::const_reference reference;
359 typedef typename vector_unary::const_pointer pointer;
361 // Construction and destruction
364 container_const_reference<self_type> (), it_ () {}
366 const_iterator (const self_type &vu, const const_subiterator_type &it):
367 container_const_reference<self_type> (vu), it_ (it) {}
371 const_iterator &operator ++ () {
376 const_iterator &operator -- () {
381 const_iterator &operator += (difference_type n) {
386 const_iterator &operator -= (difference_type n) {
391 difference_type operator - (const const_iterator &it) const {
392 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
398 const_reference operator * () const {
399 return functor_type::apply (*it_);
402 const_reference operator [] (difference_type n) const {
408 size_type index () const {
414 const_iterator &operator = (const const_iterator &it) {
415 container_const_reference<self_type>::assign (&it ());
422 bool operator == (const const_iterator &it) const {
423 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
424 return it_ == it.it_;
427 bool operator < (const const_iterator &it) const {
428 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
433 const_subiterator_type it_;
438 const_iterator begin () const {
442 const_iterator cbegin () const {
446 const_iterator end () const {
447 return find (size ());
450 const_iterator cend () const {
455 typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
458 const_reverse_iterator rbegin () const {
459 return const_reverse_iterator (end ());
462 const_reverse_iterator crbegin () const {
466 const_reverse_iterator rend () const {
467 return const_reverse_iterator (begin ());
470 const_reverse_iterator crend () const {
475 expression_closure_type e_;
478 template<class E, class F>
479 struct vector_unary_traits {
480 typedef vector_unary<E, F> expression_type;
482 // #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG
483 typedef expression_type result_type;
485 // typedef typename E::vector_temporary_type result_type;
489 // (- v) [i] = - v [i]
492 typename vector_unary_traits<E, scalar_negate<typename E::value_type> >::result_type
493 operator - (const vector_expression<E> &e) {
494 typedef typename vector_unary_traits<E, scalar_negate<typename E::value_type> >::expression_type expression_type;
495 return expression_type (e ());
498 // (conj v) [i] = conj (v [i])
501 typename vector_unary_traits<E, scalar_conj<typename E::value_type> >::result_type
502 conj (const vector_expression<E> &e) {
503 typedef typename vector_unary_traits<E, scalar_conj<typename E::value_type> >::expression_type expression_type;
504 return expression_type (e ());
507 // (real v) [i] = real (v [i])
510 typename vector_unary_traits<E, scalar_real<typename E::value_type> >::result_type
511 real (const vector_expression<E> &e) {
512 typedef typename vector_unary_traits<E, scalar_real<typename E::value_type> >::expression_type expression_type;
513 return expression_type (e ());
516 // (imag v) [i] = imag (v [i])
519 typename vector_unary_traits<E, scalar_imag<typename E::value_type> >::result_type
520 imag (const vector_expression<E> &e) {
521 typedef typename vector_unary_traits<E, scalar_imag<typename E::value_type> >::expression_type expression_type;
522 return expression_type (e ());
525 // (trans v) [i] = v [i]
528 typename vector_unary_traits<const E, scalar_identity<typename E::value_type> >::result_type
529 trans (const vector_expression<E> &e) {
530 typedef typename vector_unary_traits<const E, scalar_identity<typename E::value_type> >::expression_type expression_type;
531 return expression_type (e ());
535 typename vector_unary_traits<E, scalar_identity<typename E::value_type> >::result_type
536 trans (vector_expression<E> &e) {
537 typedef typename vector_unary_traits<E, scalar_identity<typename E::value_type> >::expression_type expression_type;
538 return expression_type (e ());
541 // (herm v) [i] = conj (v [i])
544 typename vector_unary_traits<E, scalar_conj<typename E::value_type> >::result_type
545 herm (const vector_expression<E> &e) {
546 typedef typename vector_unary_traits<E, scalar_conj<typename E::value_type> >::expression_type expression_type;
547 return expression_type (e ());
550 template<class E1, class E2, class F>
552 public vector_expression<vector_binary<E1, E2, F> > {
554 typedef E1 expression1_type;
555 typedef E2 expression2_type;
556 typedef F functor_type;
557 typedef typename E1::const_closure_type expression1_closure_type;
558 typedef typename E2::const_closure_type expression2_closure_type;
559 typedef vector_binary<E1, E2, F> self_type;
561 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
562 using vector_expression<vector_binary<E1, E2, F> >::operator ();
564 typedef typename promote_traits<typename E1::size_type, typename E2::size_type>::promote_type size_type;
565 typedef typename promote_traits<typename E1::difference_type, typename E2::difference_type>::promote_type difference_type;
566 typedef typename F::result_type value_type;
567 typedef value_type const_reference;
568 typedef const_reference reference;
569 typedef const self_type const_closure_type;
570 typedef const_closure_type closure_type;
571 typedef unknown_storage_tag storage_category;
573 // Construction and destruction
575 vector_binary (const expression1_type &e1, const expression2_type &e2):
576 e1_ (e1), e2_ (e2) {}
580 size_type size () const {
581 return BOOST_UBLAS_SAME (e1_.size (), e2_.size ());
587 const expression1_closure_type &expression1 () const {
591 const expression2_closure_type &expression2 () const {
598 const_reference operator () (size_type i) const {
599 return functor_type::apply (e1_ (i), e2_ (i));
603 const_reference operator [] (size_type i) const {
604 return functor_type::apply (e1_ [i], e2_ [i]);
607 // Closure comparison
609 bool same_closure (const vector_binary &vb) const {
610 return (*this).expression1 ().same_closure (vb.expression1 ()) &&
611 (*this).expression2 ().same_closure (vb.expression2 ());
616 typedef typename E1::const_iterator const_subiterator1_type;
617 typedef typename E2::const_iterator const_subiterator2_type;
618 typedef const value_type *const_pointer;
621 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
622 typedef typename iterator_restrict_traits<typename const_subiterator1_type::iterator_category,
623 typename const_subiterator2_type::iterator_category>::iterator_category iterator_category;
624 typedef indexed_const_iterator<const_closure_type, iterator_category> const_iterator;
625 typedef const_iterator iterator;
627 class const_iterator;
628 typedef const_iterator iterator;
633 const_iterator find (size_type i) const {
634 const_subiterator1_type it1 (e1_.find (i));
635 const_subiterator1_type it1_end (e1_.find (size ()));
636 const_subiterator2_type it2 (e2_.find (i));
637 const_subiterator2_type it2_end (e2_.find (size ()));
638 i = (std::min) (it1 != it1_end ? it1.index () : size (),
639 it2 != it2_end ? it2.index () : size ());
640 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
641 return const_iterator (*this, i);
643 return const_iterator (*this, i, it1, it1_end, it2, it2_end);
647 // Iterator merges the iterators of the referenced expressions and
648 // enhances them with the binary functor.
650 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
651 class const_iterator:
652 public container_const_reference<vector_binary>,
653 public iterator_base_traits<typename iterator_restrict_traits<typename E1::const_iterator::iterator_category,
654 typename E2::const_iterator::iterator_category>::iterator_category>::template
655 iterator_base<const_iterator, value_type>::type {
657 typedef typename iterator_restrict_traits<typename E1::const_iterator::iterator_category,
658 typename E2::const_iterator::iterator_category>::iterator_category iterator_category;
659 typedef typename vector_binary::difference_type difference_type;
660 typedef typename vector_binary::value_type value_type;
661 typedef typename vector_binary::const_reference reference;
662 typedef typename vector_binary::const_pointer pointer;
664 // Construction and destruction
667 container_const_reference<self_type> (), i_ (), it1_ (), it1_end_ (), it2_ (), it2_end_ () {}
669 const_iterator (const self_type &vb, size_type i,
670 const const_subiterator1_type &it1, const const_subiterator1_type &it1_end,
671 const const_subiterator2_type &it2, const const_subiterator2_type &it2_end):
672 container_const_reference<self_type> (vb), i_ (i), it1_ (it1), it1_end_ (it1_end), it2_ (it2), it2_end_ (it2_end) {}
675 // Dense specializations
677 void increment (dense_random_access_iterator_tag) {
678 ++ i_; ++ it1_; ++ it2_;
681 void decrement (dense_random_access_iterator_tag) {
682 -- i_; -- it1_; -- it2_;
685 void increment (dense_random_access_iterator_tag, difference_type n) {
686 i_ += n; it1_ += n; it2_ += n;
689 void decrement (dense_random_access_iterator_tag, difference_type n) {
690 i_ -= n; it1_ -= n; it2_ -= n;
693 value_type dereference (dense_random_access_iterator_tag) const {
694 return functor_type::apply (*it1_, *it2_);
697 // Packed specializations
699 void increment (packed_random_access_iterator_tag) {
700 if (it1_ != it1_end_)
701 if (it1_.index () <= i_)
703 if (it2_ != it2_end_)
704 if (it2_.index () <= i_)
709 void decrement (packed_random_access_iterator_tag) {
710 if (it1_ != it1_end_)
711 if (i_ <= it1_.index ())
713 if (it2_ != it2_end_)
714 if (i_ <= it2_.index ())
719 void increment (packed_random_access_iterator_tag, difference_type n) {
721 increment (packed_random_access_iterator_tag ());
725 decrement (packed_random_access_iterator_tag ());
730 void decrement (packed_random_access_iterator_tag, difference_type n) {
732 decrement (packed_random_access_iterator_tag ());
736 increment (packed_random_access_iterator_tag ());
741 value_type dereference (packed_random_access_iterator_tag) const {
742 typename E1::value_type t1 = typename E1::value_type/*zero*/();
743 if (it1_ != it1_end_)
744 if (it1_.index () == i_)
746 typename E2::value_type t2 = typename E2::value_type/*zero*/();
747 if (it2_ != it2_end_)
748 if (it2_.index () == i_)
750 return functor_type::apply (t1, t2);
753 // Sparse specializations
755 void increment (sparse_bidirectional_iterator_tag) {
756 size_type index1 = (*this) ().size ();
757 if (it1_ != it1_end_) {
758 if (it1_.index () <= i_)
760 if (it1_ != it1_end_)
761 index1 = it1_.index ();
763 size_type index2 = (*this) ().size ();
764 if (it2_ != it2_end_) {
765 if (it2_.index () <= i_)
767 if (it2_ != it2_end_)
768 index2 = it2_.index ();
770 i_ = (std::min) (index1, index2);
773 void decrement (sparse_bidirectional_iterator_tag) {
774 size_type index1 = (*this) ().size ();
775 if (it1_ != it1_end_) {
776 if (i_ <= it1_.index ())
778 if (it1_ != it1_end_)
779 index1 = it1_.index ();
781 size_type index2 = (*this) ().size ();
782 if (it2_ != it2_end_) {
783 if (i_ <= it2_.index ())
785 if (it2_ != it2_end_)
786 index2 = it2_.index ();
788 i_ = (std::max) (index1, index2);
791 void increment (sparse_bidirectional_iterator_tag, difference_type n) {
793 increment (sparse_bidirectional_iterator_tag ());
797 decrement (sparse_bidirectional_iterator_tag ());
802 void decrement (sparse_bidirectional_iterator_tag, difference_type n) {
804 decrement (sparse_bidirectional_iterator_tag ());
808 increment (sparse_bidirectional_iterator_tag ());
813 value_type dereference (sparse_bidirectional_iterator_tag) const {
814 typename E1::value_type t1 = typename E1::value_type/*zero*/();
815 if (it1_ != it1_end_)
816 if (it1_.index () == i_)
818 typename E2::value_type t2 = typename E2::value_type/*zero*/();
819 if (it2_ != it2_end_)
820 if (it2_.index () == i_)
822 return static_cast<value_type>(functor_type::apply (t1, t2));
828 const_iterator &operator ++ () {
829 increment (iterator_category ());
833 const_iterator &operator -- () {
834 decrement (iterator_category ());
838 const_iterator &operator += (difference_type n) {
839 increment (iterator_category (), n);
843 const_iterator &operator -= (difference_type n) {
844 decrement (iterator_category (), n);
848 difference_type operator - (const const_iterator &it) const {
849 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
850 return index () - it.index ();
855 const_reference operator * () const {
856 return dereference (iterator_category ());
859 const_reference operator [] (difference_type n) const {
865 size_type index () const {
871 const_iterator &operator = (const const_iterator &it) {
872 container_const_reference<self_type>::assign (&it ());
875 it1_end_ = it.it1_end_;
877 it2_end_ = it.it2_end_;
883 bool operator == (const const_iterator &it) const {
884 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
885 return index () == it.index ();
888 bool operator < (const const_iterator &it) const {
889 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
890 return index () < it.index ();
895 const_subiterator1_type it1_;
896 const_subiterator1_type it1_end_;
897 const_subiterator2_type it2_;
898 const_subiterator2_type it2_end_;
903 const_iterator begin () const {
907 const_iterator cbegin () const {
911 const_iterator end () const {
912 return find (size ());
915 const_iterator cend () const {
920 typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
923 const_reverse_iterator rbegin () const {
924 return const_reverse_iterator (end ());
927 const_reverse_iterator crbegin () const {
931 const_reverse_iterator rend () const {
932 return const_reverse_iterator (begin ());
935 const_reverse_iterator crend () const {
940 expression1_closure_type e1_;
941 expression2_closure_type e2_;
944 template<class E1, class E2, class F>
945 struct vector_binary_traits {
946 typedef vector_binary<E1, E2, F> expression_type;
947 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG
948 typedef expression_type result_type;
950 typedef typename E1::vector_temporary_type result_type;
954 // (v1 + v2) [i] = v1 [i] + v2 [i]
955 template<class E1, class E2>
957 typename vector_binary_traits<E1, E2, scalar_plus<typename E1::value_type,
958 typename E2::value_type> >::result_type
959 operator + (const vector_expression<E1> &e1,
960 const vector_expression<E2> &e2) {
961 typedef typename vector_binary_traits<E1, E2, scalar_plus<typename E1::value_type,
962 typename E2::value_type> >::expression_type expression_type;
963 return expression_type (e1 (), e2 ());
966 // (v1 - v2) [i] = v1 [i] - v2 [i]
967 template<class E1, class E2>
969 typename vector_binary_traits<E1, E2, scalar_minus<typename E1::value_type,
970 typename E2::value_type> >::result_type
971 operator - (const vector_expression<E1> &e1,
972 const vector_expression<E2> &e2) {
973 typedef typename vector_binary_traits<E1, E2, scalar_minus<typename E1::value_type,
974 typename E2::value_type> >::expression_type expression_type;
975 return expression_type (e1 (), e2 ());
978 // (v1 * v2) [i] = v1 [i] * v2 [i]
979 template<class E1, class E2>
981 typename vector_binary_traits<E1, E2, scalar_multiplies<typename E1::value_type,
982 typename E2::value_type> >::result_type
983 element_prod (const vector_expression<E1> &e1,
984 const vector_expression<E2> &e2) {
985 typedef typename vector_binary_traits<E1, E2, scalar_multiplies<typename E1::value_type,
986 typename E2::value_type> >::expression_type expression_type;
987 return expression_type (e1 (), e2 ());
990 // (v1 / v2) [i] = v1 [i] / v2 [i]
991 template<class E1, class E2>
993 typename vector_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
994 typename E2::value_type> >::result_type
995 element_div (const vector_expression<E1> &e1,
996 const vector_expression<E2> &e2) {
997 typedef typename vector_binary_traits<E1, E2, scalar_divides<typename E1::value_type,
998 typename E2::value_type> >::expression_type expression_type;
999 return expression_type (e1 (), e2 ());
1003 template<class E1, class E2, class F>
1004 class vector_binary_scalar1:
1005 public vector_expression<vector_binary_scalar1<E1, E2, F> > {
1007 typedef F functor_type;
1008 typedef E1 expression1_type;
1009 typedef E2 expression2_type;
1011 typedef const E1& expression1_closure_type;
1012 typedef typename E2::const_closure_type expression2_closure_type;
1014 typedef vector_binary_scalar1<E1, E2, F> self_type;
1016 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1017 using vector_expression<vector_binary_scalar1<E1, E2, F> >::operator ();
1019 typedef typename E2::size_type size_type;
1020 typedef typename E2::difference_type difference_type;
1021 typedef typename F::result_type value_type;
1022 typedef value_type const_reference;
1023 typedef const_reference reference;
1024 typedef const self_type const_closure_type;
1025 typedef const_closure_type closure_type;
1026 typedef unknown_storage_tag storage_category;
1028 // Construction and destruction
1030 vector_binary_scalar1 (const expression1_type &e1, const expression2_type &e2):
1031 e1_ (e1), e2_ (e2) {}
1035 size_type size () const {
1042 const_reference operator () (size_type i) const {
1043 return functor_type::apply (e1_, e2_ (i));
1047 const_reference operator [] (size_type i) const {
1048 return functor_type::apply (e1_, e2_ [i]);
1051 // Closure comparison
1053 bool same_closure (const vector_binary_scalar1 &vbs1) const {
1054 return &e1_ == &(vbs1.e1_) &&
1055 (*this).e2_.same_closure (vbs1.e2_);
1060 typedef expression1_type const_subiterator1_type;
1061 typedef typename expression2_type::const_iterator const_subiterator2_type;
1062 typedef const value_type *const_pointer;
1065 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1066 typedef indexed_const_iterator<const_closure_type, typename const_subiterator2_type::iterator_category> const_iterator;
1067 typedef const_iterator iterator;
1069 class const_iterator;
1070 typedef const_iterator iterator;
1075 const_iterator find (size_type i) const {
1076 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1077 const_subiterator2_type it (e2_.find (i));
1078 return const_iterator (*this, it.index ());
1080 return const_iterator (*this, const_subiterator1_type (e1_), e2_.find (i));
1084 // Iterator enhances the iterator of the referenced vector expression
1085 // with the binary functor.
1087 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1088 class const_iterator:
1089 public container_const_reference<vector_binary_scalar1>,
1090 public iterator_base_traits<typename E2::const_iterator::iterator_category>::template
1091 iterator_base<const_iterator, value_type>::type {
1093 typedef typename E2::const_iterator::iterator_category iterator_category;
1094 typedef typename vector_binary_scalar1::difference_type difference_type;
1095 typedef typename vector_binary_scalar1::value_type value_type;
1096 typedef typename vector_binary_scalar1::const_reference reference;
1097 typedef typename vector_binary_scalar1::const_pointer pointer;
1099 // Construction and destruction
1102 container_const_reference<self_type> (), it1_ (), it2_ () {}
1104 const_iterator (const self_type &vbs, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
1105 container_const_reference<self_type> (vbs), it1_ (it1), it2_ (it2) {}
1109 const_iterator &operator ++ () {
1114 const_iterator &operator -- () {
1119 const_iterator &operator += (difference_type n) {
1124 const_iterator &operator -= (difference_type n) {
1129 difference_type operator - (const const_iterator &it) const {
1130 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1131 // FIXME we shouldn't compare floats
1132 // BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1133 return it2_ - it.it2_;
1138 const_reference operator * () const {
1139 return functor_type::apply (it1_, *it2_);
1142 const_reference operator [] (difference_type n) const {
1143 return *(*this + n);
1148 size_type index () const {
1149 return it2_.index ();
1154 const_iterator &operator = (const const_iterator &it) {
1155 container_const_reference<self_type>::assign (&it ());
1163 bool operator == (const const_iterator &it) const {
1164 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1165 // FIXME we shouldn't compare floats
1166 // BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1167 return it2_ == it.it2_;
1170 bool operator < (const const_iterator &it) const {
1171 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1172 // FIXME we shouldn't compare floats
1173 // BOOST_UBLAS_CHECK (it1_ == it.it1_, external_logic ());
1174 return it2_ < it.it2_;
1178 const_subiterator1_type it1_;
1179 const_subiterator2_type it2_;
1184 const_iterator begin () const {
1188 const_iterator cbegin () const {
1192 const_iterator end () const {
1193 return find (size ());
1196 const_iterator cend () const {
1201 typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
1204 const_reverse_iterator rbegin () const {
1205 return const_reverse_iterator (end ());
1208 const_reverse_iterator crbegin () const {
1212 const_reverse_iterator rend () const {
1213 return const_reverse_iterator (begin ());
1216 const_reverse_iterator crend () const {
1221 expression1_closure_type e1_;
1222 expression2_closure_type e2_;
1225 template<class E1, class E2, class F>
1226 struct vector_binary_scalar1_traits {
1227 typedef vector_binary_scalar1<E1, E2, F> expression_type; // allow E1 to be builtin type
1228 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG
1229 typedef expression_type result_type;
1231 typedef typename E2::vector_temporary_type result_type;
1235 // (t * v) [i] = t * v [i]
1236 template<class T1, class E2>
1238 typename boost::enable_if< is_convertible<T1, typename E2::value_type >,
1239 typename vector_binary_scalar1_traits<const T1, E2, scalar_multiplies<T1, typename E2::value_type> >::result_type
1241 operator * (const T1 &e1,
1242 const vector_expression<E2> &e2) {
1243 typedef typename vector_binary_scalar1_traits<const T1, E2, scalar_multiplies<T1, typename E2::value_type> >::expression_type expression_type;
1244 return expression_type (e1, e2 ());
1248 template<class E1, class E2, class F>
1249 class vector_binary_scalar2:
1250 public vector_expression<vector_binary_scalar2<E1, E2, F> > {
1252 typedef F functor_type;
1253 typedef E1 expression1_type;
1254 typedef E2 expression2_type;
1255 typedef typename E1::const_closure_type expression1_closure_type;
1256 typedef const E2& expression2_closure_type;
1257 typedef vector_binary_scalar2<E1, E2, F> self_type;
1259 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
1260 using vector_expression<vector_binary_scalar2<E1, E2, F> >::operator ();
1262 typedef typename E1::size_type size_type;
1263 typedef typename E1::difference_type difference_type;
1264 typedef typename F::result_type value_type;
1265 typedef value_type const_reference;
1266 typedef const_reference reference;
1267 typedef const self_type const_closure_type;
1268 typedef const_closure_type closure_type;
1269 typedef unknown_storage_tag storage_category;
1271 // Construction and destruction
1273 vector_binary_scalar2 (const expression1_type &e1, const expression2_type &e2):
1274 e1_ (e1), e2_ (e2) {}
1278 size_type size () const {
1285 const_reference operator () (size_type i) const {
1286 return functor_type::apply (e1_ (i), e2_);
1290 const_reference operator [] (size_type i) const {
1291 return functor_type::apply (e1_ [i], e2_);
1294 // Closure comparison
1296 bool same_closure (const vector_binary_scalar2 &vbs2) const {
1297 return (*this).e1_.same_closure (vbs2.e1_) &&
1298 &e2_ == &(vbs2.e2_);
1303 typedef typename expression1_type::const_iterator const_subiterator1_type;
1304 typedef expression2_type const_subiterator2_type;
1305 typedef const value_type *const_pointer;
1308 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1309 typedef indexed_const_iterator<const_closure_type, typename const_subiterator2_type::iterator_category> const_iterator;
1310 typedef const_iterator iterator;
1312 class const_iterator;
1313 typedef const_iterator iterator;
1318 const_iterator find (size_type i) const {
1319 #ifdef BOOST_UBLAS_USE_INDEXED_ITERATOR
1320 const_subiterator1_type it (e1_.find (i));
1321 return const_iterator (*this, it.index ());
1323 return const_iterator (*this, e1_.find (i), const_subiterator2_type (e2_));
1327 // Iterator enhances the iterator of the referenced vector expression
1328 // with the binary functor.
1330 #ifndef BOOST_UBLAS_USE_INDEXED_ITERATOR
1331 class const_iterator:
1332 public container_const_reference<vector_binary_scalar2>,
1333 public iterator_base_traits<typename E1::const_iterator::iterator_category>::template
1334 iterator_base<const_iterator, value_type>::type {
1336 typedef typename E1::const_iterator::iterator_category iterator_category;
1337 typedef typename vector_binary_scalar2::difference_type difference_type;
1338 typedef typename vector_binary_scalar2::value_type value_type;
1339 typedef typename vector_binary_scalar2::const_reference reference;
1340 typedef typename vector_binary_scalar2::const_pointer pointer;
1342 // Construction and destruction
1345 container_const_reference<self_type> (), it1_ (), it2_ () {}
1347 const_iterator (const self_type &vbs, const const_subiterator1_type &it1, const const_subiterator2_type &it2):
1348 container_const_reference<self_type> (vbs), it1_ (it1), it2_ (it2) {}
1352 const_iterator &operator ++ () {
1357 const_iterator &operator -- () {
1362 const_iterator &operator += (difference_type n) {
1367 const_iterator &operator -= (difference_type n) {
1372 difference_type operator - (const const_iterator &it) const {
1373 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1374 // FIXME we shouldn't compare floats
1375 // BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
1376 return it1_ - it.it1_;
1381 const_reference operator * () const {
1382 return functor_type::apply (*it1_, it2_);
1385 const_reference operator [] (difference_type n) const {
1386 return *(*this + n);
1391 size_type index () const {
1392 return it1_.index ();
1397 const_iterator &operator = (const const_iterator &it) {
1398 container_const_reference<self_type>::assign (&it ());
1406 bool operator == (const const_iterator &it) const {
1407 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1408 // FIXME we shouldn't compare floats
1409 // BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
1410 return it1_ == it.it1_;
1413 bool operator < (const const_iterator &it) const {
1414 BOOST_UBLAS_CHECK ((*this) ().same_closure (it ()), external_logic ());
1415 // FIXME we shouldn't compare floats
1416 // BOOST_UBLAS_CHECK (it2_ == it.it2_, external_logic ());
1417 return it1_ < it.it1_;
1421 const_subiterator1_type it1_;
1422 const_subiterator2_type it2_;
1427 const_iterator begin () const {
1431 const_iterator cbegin () const {
1435 const_iterator end () const {
1436 return find (size ());
1439 const_iterator cend () const {
1444 typedef reverse_iterator_base<const_iterator> const_reverse_iterator;
1447 const_reverse_iterator rbegin () const {
1448 return const_reverse_iterator (end ());
1451 const_reverse_iterator crbegin () const {
1455 const_reverse_iterator rend () const {
1456 return const_reverse_iterator (begin ());
1459 const_reverse_iterator crend () const {
1464 expression1_closure_type e1_;
1465 expression2_closure_type e2_;
1468 template<class E1, class E2, class F>
1469 struct vector_binary_scalar2_traits {
1470 typedef vector_binary_scalar2<E1, E2, F> expression_type; // allow E2 to be builtin type
1471 #ifndef BOOST_UBLAS_SIMPLE_ET_DEBUG
1472 typedef expression_type result_type;
1474 typedef typename E1::vector_temporary_type result_type;
1478 // (v * t) [i] = v [i] * t
1479 template<class E1, class T2>
1481 typename boost::enable_if< is_convertible<T2, typename E1::value_type >,
1482 typename vector_binary_scalar2_traits<E1, const T2, scalar_multiplies<typename E1::value_type, T2> >::result_type
1484 operator * (const vector_expression<E1> &e1,
1486 typedef typename vector_binary_scalar2_traits<E1, const T2, scalar_multiplies<typename E1::value_type, T2> >::expression_type expression_type;
1487 return expression_type (e1 (), e2);
1490 // (v / t) [i] = v [i] / t
1491 template<class E1, class T2>
1493 typename boost::enable_if< is_convertible<T2, typename E1::value_type >,
1494 typename vector_binary_scalar2_traits<E1, const T2, scalar_divides<typename E1::value_type, T2> >::result_type
1496 operator / (const vector_expression<E1> &e1,
1498 typedef typename vector_binary_scalar2_traits<E1, const T2, scalar_divides<typename E1::value_type, T2> >::expression_type expression_type;
1499 return expression_type (e1 (), e2);
1503 template<class E, class F>
1504 class vector_scalar_unary:
1505 public scalar_expression<vector_scalar_unary<E, F> > {
1507 typedef E expression_type;
1508 typedef F functor_type;
1509 typedef typename E::const_closure_type expression_closure_type;
1510 typedef typename E::const_iterator::iterator_category iterator_category;
1511 typedef vector_scalar_unary<E, F> self_type;
1513 typedef typename F::result_type value_type;
1514 typedef typename E::difference_type difference_type;
1515 typedef const self_type const_closure_type;
1516 typedef const_closure_type closure_type;
1517 typedef unknown_storage_tag storage_category;
1519 // Construction and destruction
1521 explicit vector_scalar_unary (const expression_type &e):
1525 // Expression accessors
1527 const expression_closure_type &expression () const {
1533 operator value_type () const {
1534 return evaluate (iterator_category ());
1538 // Dense random access specialization
1540 value_type evaluate (dense_random_access_iterator_tag) const {
1541 #ifdef BOOST_UBLAS_USE_INDEXING
1542 return functor_type::apply (e_);
1543 #elif BOOST_UBLAS_USE_ITERATING
1544 difference_type size = e_.size ();
1545 return functor_type::apply (size, e_.begin ());
1547 difference_type size = e_.size ();
1548 if (size >= BOOST_UBLAS_ITERATOR_THRESHOLD)
1549 return functor_type::apply (size, e_.begin ());
1551 return functor_type::apply (e_);
1555 // Packed bidirectional specialization
1557 value_type evaluate (packed_random_access_iterator_tag) const {
1558 return functor_type::apply (e_.begin (), e_.end ());
1561 // Sparse bidirectional specialization
1563 value_type evaluate (sparse_bidirectional_iterator_tag) const {
1564 return functor_type::apply (e_.begin (), e_.end ());
1568 expression_closure_type e_;
1571 template<class E, class F>
1572 struct vector_scalar_unary_traits {
1573 typedef vector_scalar_unary<E, F> expression_type;
1574 #if !defined (BOOST_UBLAS_SIMPLE_ET_DEBUG) && defined (BOOST_UBLAS_USE_SCALAR_ET)
1575 // FIXME don't define USE_SCALAR_ET other then for testing
1576 // They do not work for complex types
1577 typedef expression_type result_type;
1579 typedef typename F::result_type result_type;
1583 // sum v = sum (v [i])
1586 typename vector_scalar_unary_traits<E, vector_sum<E> >::result_type
1587 sum (const vector_expression<E> &e) {
1588 typedef typename vector_scalar_unary_traits<E, vector_sum<E> >::expression_type expression_type;
1589 return expression_type (e ());
1592 // real: norm_1 v = sum (abs (v [i]))
1593 // complex: norm_1 v = sum (abs (real (v [i])) + abs (imag (v [i])))
1596 typename vector_scalar_unary_traits<E, vector_norm_1<E> >::result_type
1597 norm_1 (const vector_expression<E> &e) {
1598 typedef typename vector_scalar_unary_traits<E, vector_norm_1<E> >::expression_type expression_type;
1599 return expression_type (e ());
1602 // real: norm_2 v = sqrt (sum (v [i] * v [i]))
1603 // complex: norm_2 v = sqrt (sum (v [i] * conj (v [i])))
1606 typename vector_scalar_unary_traits<E, vector_norm_2<E> >::result_type
1607 norm_2 (const vector_expression<E> &e) {
1608 typedef typename vector_scalar_unary_traits<E, vector_norm_2<E> >::expression_type expression_type;
1609 return expression_type (e ());
1612 // real: norm_2_square v = sum(v [i] * v [i])
1613 // complex: norm_2_square v = sum(v [i] * conj (v [i]))
1616 typename vector_scalar_unary_traits<E, vector_norm_2_square<E> >::result_type
1617 norm_2_square (const vector_expression<E> &e) {
1618 typedef typename vector_scalar_unary_traits<E, vector_norm_2_square<E> >::expression_type expression_type;
1619 return expression_type (e ());
1622 // real: norm_inf v = maximum (abs (v [i]))
1623 // complex: norm_inf v = maximum (maximum (abs (real (v [i])), abs (imag (v [i]))))
1626 typename vector_scalar_unary_traits<E, vector_norm_inf<E> >::result_type
1627 norm_inf (const vector_expression<E> &e) {
1628 typedef typename vector_scalar_unary_traits<E, vector_norm_inf<E> >::expression_type expression_type;
1629 return expression_type (e ());
1632 // real: index_norm_inf v = minimum (i: abs (v [i]) == maximum (abs (v [i])))
1635 typename vector_scalar_unary_traits<E, vector_index_norm_inf<E> >::result_type
1636 index_norm_inf (const vector_expression<E> &e) {
1637 typedef typename vector_scalar_unary_traits<E, vector_index_norm_inf<E> >::expression_type expression_type;
1638 return expression_type (e ());
1641 template<class E1, class E2, class F>
1642 class vector_scalar_binary:
1643 public scalar_expression<vector_scalar_binary<E1, E2, F> > {
1645 typedef E1 expression1_type;
1646 typedef E2 expression2_type;
1647 typedef F functor_type;
1648 typedef typename E1::const_closure_type expression1_closure_type;
1649 typedef typename E2::const_closure_type expression2_closure_type;
1650 typedef typename iterator_restrict_traits<typename E1::const_iterator::iterator_category,
1651 typename E2::const_iterator::iterator_category>::iterator_category iterator_category;
1652 typedef vector_scalar_binary<E1, E2, F> self_type;
1654 static const unsigned complexity = 1;
1655 typedef typename F::result_type value_type;
1656 typedef typename E1::difference_type difference_type;
1657 typedef const self_type const_closure_type;
1658 typedef const_closure_type closure_type;
1659 typedef unknown_storage_tag storage_category;
1661 // Construction and destruction
1663 vector_scalar_binary (const expression1_type &e1, const expression2_type &e2):
1664 e1_ (e1), e2_ (e2) {}
1669 const expression1_closure_type &expression1 () const {
1673 const expression2_closure_type &expression2 () const {
1679 operator value_type () const {
1680 return evaluate (iterator_category ());
1684 // Dense random access specialization
1686 value_type evaluate (dense_random_access_iterator_tag) const {
1687 BOOST_UBLAS_CHECK (e1_.size () == e2_.size (), external_logic());
1688 #ifdef BOOST_UBLAS_USE_INDEXING
1689 return functor_type::apply (e1_, e2_);
1690 #elif BOOST_UBLAS_USE_ITERATING
1691 difference_type size = BOOST_UBLAS_SAME (e1_.size (), e2_.size ());
1692 return functor_type::apply (size, e1_.begin (), e2_.begin ());
1694 difference_type size = BOOST_UBLAS_SAME (e1_.size (), e2_.size ());
1695 if (size >= BOOST_UBLAS_ITERATOR_THRESHOLD)
1696 return functor_type::apply (size, e1_.begin (), e2_.begin ());
1698 return functor_type::apply (e1_, e2_);
1702 // Packed bidirectional specialization
1704 value_type evaluate (packed_random_access_iterator_tag) const {
1705 BOOST_UBLAS_CHECK (e1_.size () == e2_.size (), external_logic());
1706 return functor_type::apply (e1_.begin (), e1_.end (), e2_.begin (), e2_.end ());
1709 // Sparse bidirectional specialization
1711 value_type evaluate (sparse_bidirectional_iterator_tag) const {
1712 BOOST_UBLAS_CHECK (e1_.size () == e2_.size (), external_logic());
1713 return functor_type::apply (e1_.begin (), e1_.end (), e2_.begin (), e2_.end (), sparse_bidirectional_iterator_tag ());
1717 expression1_closure_type e1_;
1718 expression2_closure_type e2_;
1721 template<class E1, class E2, class F>
1722 struct vector_scalar_binary_traits {
1723 typedef vector_scalar_binary<E1, E2, F> expression_type;
1724 #if !defined (BOOST_UBLAS_SIMPLE_ET_DEBUG) && defined (BOOST_UBLAS_USE_SCALAR_ET)
1725 // FIXME don't define USE_SCALAR_ET other then for testing
1726 // They do not work for complex types
1727 typedef expression_type result_type;
1729 typedef typename F::result_type result_type;
1733 // inner_prod (v1, v2) = sum (v1 [i] * v2 [i])
1734 template<class E1, class E2>
1736 typename vector_scalar_binary_traits<E1, E2, vector_inner_prod<E1, E2,
1737 typename promote_traits<typename E1::value_type,
1738 typename E2::value_type>::promote_type> >::result_type
1739 inner_prod (const vector_expression<E1> &e1,
1740 const vector_expression<E2> &e2) {
1741 typedef typename vector_scalar_binary_traits<E1, E2, vector_inner_prod<E1, E2,
1742 typename promote_traits<typename E1::value_type,
1743 typename E2::value_type>::promote_type> >::expression_type expression_type;
1744 return expression_type (e1 (), e2 ());
1747 template<class E1, class E2>
1749 typename vector_scalar_binary_traits<E1, E2, vector_inner_prod<E1, E2,
1750 typename type_traits<typename promote_traits<typename E1::value_type,
1751 typename E2::value_type>::promote_type>::precision_type> >::result_type
1752 prec_inner_prod (const vector_expression<E1> &e1,
1753 const vector_expression<E2> &e2) {
1754 typedef typename vector_scalar_binary_traits<E1, E2, vector_inner_prod<E1, E2,
1755 typename type_traits<typename promote_traits<typename E1::value_type,
1756 typename E2::value_type>::promote_type>::precision_type> >::expression_type expression_type;
1757 return expression_type (e1 (), e2 ());