Boost.uBlas 1.49
Linear Algebra in C++: matrices, vectors and numeric algorithms

sparse_view.hpp

Go to the documentation of this file.
00001 //
00002 //  Copyright (c) 2009
00003 //  Gunter Winkler
00004 //
00005 //  Distributed under the Boost Software License, Version 1.0. (See
00006 //  accompanying file LICENSE_1_0.txt or copy at
00007 //  http://www.boost.org/LICENSE_1_0.txt)
00008 //
00009 //
00010 
00011 #ifndef _BOOST_UBLAS_SPARSE_VIEW_
00012 #define _BOOST_UBLAS_SPARSE_VIEW_
00013 
00014 #include <boost/numeric/ublas/matrix_expression.hpp>
00015 #include <boost/numeric/ublas/detail/matrix_assign.hpp>
00016 #if BOOST_UBLAS_TYPE_CHECK
00017 #include <boost/numeric/ublas/matrix.hpp>
00018 #endif
00019 
00020 #include <boost/next_prior.hpp>
00021 #include <boost/type_traits/remove_cv.hpp>
00022 
00023 namespace boost { namespace numeric { namespace ublas {
00024 
00025     // view a chunk of memory as ublas array
00026 
00027     template < class T >
00028     class c_array_view
00029         : public storage_array< c_array_view<T> > {
00030     private:
00031         typedef c_array_view<T> self_type;
00032         typedef T * pointer;
00033 
00034     public:
00035         // TODO: think about a const pointer
00036         typedef const pointer array_type;
00037 
00038         typedef std::size_t size_type;
00039         typedef std::ptrdiff_t difference_type;
00040 
00041         typedef T value_type;
00042         typedef const T  &const_reference;
00043         typedef const T  *const_pointer;
00044 
00045         typedef const_pointer const_iterator;
00046         typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
00047 
00048         //
00049         // typedefs required by vector concept
00050         //
00051 
00052         typedef dense_tag  storage_category;
00053         typedef const vector_reference<const self_type>    const_closure_type;
00054 
00055         c_array_view(size_type size, array_type data) :
00056             size_(size), data_(data)
00057         {}
00058 
00059         ~c_array_view()
00060         {}
00061 
00062         //
00063         // immutable methods of container concept
00064         //
00065 
00066         BOOST_UBLAS_INLINE
00067         size_type size () const {
00068             return size_;
00069         }
00070 
00071         BOOST_UBLAS_INLINE
00072         const_reference operator [] (size_type i) const {
00073             BOOST_UBLAS_CHECK (i < size_, bad_index ());
00074             return data_ [i];
00075         }
00076 
00077         BOOST_UBLAS_INLINE
00078         const_iterator begin () const {
00079             return data_;
00080         }
00081         BOOST_UBLAS_INLINE
00082         const_iterator end () const {
00083             return data_ + size_;
00084         }
00085 
00086         BOOST_UBLAS_INLINE
00087         const_reverse_iterator rbegin () const {
00088             return const_reverse_iterator (end ());
00089         }
00090         BOOST_UBLAS_INLINE
00091         const_reverse_iterator rend () const {
00092             return const_reverse_iterator (begin ());
00093         }
00094 
00095     private:
00096         size_type  size_;
00097         array_type data_;
00098     };
00099 
00100 
00115     template<class L, std::size_t IB, class IA, class JA, class TA>
00116     class compressed_matrix_view:
00117         public matrix_expression<compressed_matrix_view<L, IB, IA, JA, TA> > {
00118 
00119     public:
00120         typedef typename vector_view_traits<TA>::value_type value_type;
00121 
00122     private:
00123         typedef value_type &true_reference;
00124         typedef value_type *pointer;
00125         typedef const value_type *const_pointer;
00126         typedef L layout_type;
00127         typedef compressed_matrix_view<L, IB, IA, JA, TA> self_type;
00128 
00129     public:
00130 #ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
00131         using matrix_expression<self_type>::operator ();
00132 #endif
00133         // ISSUE require type consistency check
00134         // is_convertable (IA::size_type, TA::size_type)
00135         typedef typename boost::remove_cv<typename vector_view_traits<JA>::value_type>::type index_type;
00136         // for compatibility, should be removed some day ...
00137         typedef index_type size_type;
00138         // size_type for the data arrays.
00139         typedef typename vector_view_traits<JA>::size_type array_size_type;
00140         typedef typename vector_view_traits<JA>::difference_type difference_type;
00141         typedef const value_type & const_reference;
00142 
00143         // do NOT define reference type, because class is read only
00144         // typedef value_type & reference;
00145 
00146         typedef IA rowptr_array_type;
00147         typedef JA index_array_type;
00148         typedef TA value_array_type;
00149         typedef const matrix_reference<const self_type> const_closure_type;
00150         typedef matrix_reference<self_type> closure_type;
00151 
00152         // FIXME: define a corresponding temporary type
00153         // typedef compressed_vector<T, IB, IA, TA> vector_temporary_type;
00154 
00155         // FIXME: define a corresponding temporary type
00156         // typedef self_type matrix_temporary_type;
00157 
00158         typedef sparse_tag storage_category;
00159         typedef typename L::orientation_category orientation_category;
00160 
00161         //
00162         // private types for internal use
00163         //
00164 
00165     private:
00166         typedef typename vector_view_traits<index_array_type>::const_iterator const_subiterator_type;
00167 
00168         //
00169         // Construction and destruction
00170         //
00171     private:
00173         BOOST_UBLAS_INLINE
00174         compressed_matrix_view () { }
00175 
00176     public:
00177         BOOST_UBLAS_INLINE
00178         compressed_matrix_view (index_type n_rows, index_type n_cols, array_size_type nnz
00179                                 , const rowptr_array_type & iptr
00180                                 , const index_array_type & jptr
00181                                 , const value_array_type & values):
00182             matrix_expression<self_type> (),
00183             size1_ (n_rows), size2_ (n_cols), 
00184             nnz_ (nnz),
00185             index1_data_ (iptr), 
00186             index2_data_ (jptr), 
00187             value_data_ (values) {
00188             storage_invariants ();
00189         }
00190 
00191         BOOST_UBLAS_INLINE
00192         compressed_matrix_view(const compressed_matrix_view& o) :
00193             size1_(size1_), size2_(size2_),
00194             nnz_(nnz_),
00195             index1_data_(index1_data_),
00196             index2_data_(index2_data_),
00197             value_data_(value_data_)
00198         {}
00199 
00200         //
00201         // implement immutable iterator types
00202         //
00203 
00204         class const_iterator1 {};
00205         class const_iterator2 {};
00206 
00207         typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
00208         typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
00209 
00210         //
00211         // implement all read only methods for the matrix expression concept
00212         // 
00213 
00215         index_type size1() const {
00216             return size1_;
00217         }
00218 
00220         index_type size2() const {
00221             return size2_;
00222         }
00223 
00225         value_type operator()(index_type i, index_type j) const {
00226             const_pointer p = find_element(i,j);
00227             if (!p) {
00228                 return zero_;
00229             } else {
00230                 return *p;
00231             }
00232         }
00233         
00234 
00235     private:
00236         //
00237         // private helper functions
00238         //
00239 
00240         const_pointer find_element (index_type i, index_type j) const {
00241             index_type element1 (layout_type::index_M (i, j));
00242             index_type element2 (layout_type::index_m (i, j));
00243 
00244             const array_size_type itv      = zero_based( index1_data_[element1] );
00245             const array_size_type itv_next = zero_based( index1_data_[element1+1] );
00246 
00247             const_subiterator_type it_start = boost::next(vector_view_traits<index_array_type>::begin(index2_data_),itv);
00248             const_subiterator_type it_end = boost::next(vector_view_traits<index_array_type>::begin(index2_data_),itv_next);
00249             const_subiterator_type it = find_index_in_row(it_start, it_end, element2) ;
00250             
00251             if (it == it_end || *it != k_based (element2))
00252                 return 0;
00253             return &value_data_ [it - vector_view_traits<index_array_type>::begin(index2_data_)];
00254         }
00255 
00256         const_subiterator_type find_index_in_row(const_subiterator_type it_start
00257                                                  , const_subiterator_type it_end
00258                                                  , index_type index) const {
00259             return std::lower_bound( it_start
00260                                      , it_end
00261                                      , k_based (index) );
00262         }
00263 
00264 
00265     private:
00266         void storage_invariants () const {
00267             BOOST_UBLAS_CHECK (index1_data_ [layout_type::size_M (size1_, size2_)] == k_based (nnz_), external_logic ());
00268         }
00269         
00270         index_type size1_;
00271         index_type size2_;
00272 
00273         array_size_type nnz_;
00274 
00275         const rowptr_array_type & index1_data_;
00276         const index_array_type & index2_data_;
00277         const value_array_type & value_data_;
00278 
00279         static const value_type zero_;
00280 
00281         BOOST_UBLAS_INLINE
00282         static index_type zero_based (index_type k_based_index) {
00283             return k_based_index - IB;
00284         }
00285         BOOST_UBLAS_INLINE
00286         static index_type k_based (index_type zero_based_index) {
00287             return zero_based_index + IB;
00288         }
00289 
00290         friend class iterator1;
00291         friend class iterator2;
00292         friend class const_iterator1;
00293         friend class const_iterator2;
00294     };
00295 
00296     template<class L, std::size_t IB, class IA, class JA, class TA  >
00297     const typename compressed_matrix_view<L,IB,IA,JA,TA>::value_type 
00298     compressed_matrix_view<L,IB,IA,JA,TA>::zero_ = value_type/*zero*/();
00299 
00300 
00301     template<class L, std::size_t IB, class IA, class JA, class TA  >
00302     compressed_matrix_view<L,IB,IA,JA,TA>
00303     make_compressed_matrix_view(typename vector_view_traits<JA>::value_type n_rows
00304                                 , typename vector_view_traits<JA>::value_type n_cols
00305                                 , typename vector_view_traits<JA>::size_type nnz
00306                                 , const IA & ia
00307                                 , const JA & ja
00308                                 , const TA & ta) {
00309 
00310         return compressed_matrix_view<L,IB,IA,JA,TA>(n_rows, n_cols, nnz, ia, ja, ta);
00311 
00312     }
00313 
00314 }}}
00315 
00316 #endif