Planet
navi homePPSaboutscreenshotsdownloaddevelopmentforum

source: downloads/boost_1_34_1/boost/numeric/ublas/vector_expression.hpp @ 33

Last change on this file since 33 was 29, checked in by landauf, 16 years ago

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