![]() |
Boost.uBlas 1.49
Linear Algebra in C++: matrices, vectors and numeric algorithms
|
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