New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
sparse_view.hpp in vendors/XIOS/current/extern/boost/include/boost/numeric/ublas/experimental – NEMO

source: vendors/XIOS/current/extern/boost/include/boost/numeric/ublas/experimental/sparse_view.hpp @ 3428

Last change on this file since 3428 was 3428, checked in by rblod, 12 years ago

importing initial XIOS vendor drop

File size: 10.2 KB
Line 
1//
2//  Copyright (c) 2009
3//  Gunter Winkler
4//
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)
8//
9//
10
11#ifndef _BOOST_UBLAS_SPARSE_VIEW_
12#define _BOOST_UBLAS_SPARSE_VIEW_
13
14#include <boost/numeric/ublas/matrix_expression.hpp>
15#include <boost/numeric/ublas/detail/matrix_assign.hpp>
16#if BOOST_UBLAS_TYPE_CHECK
17#include <boost/numeric/ublas/matrix.hpp>
18#endif
19
20#include <boost/next_prior.hpp>
21#include <boost/type_traits/remove_cv.hpp>
22
23namespace boost { namespace numeric { namespace ublas {
24
25    // view a chunk of memory as ublas array
26
27    template < class T >
28    class c_array_view
29        : public storage_array< c_array_view<T> > {
30    private:
31        typedef c_array_view<T> self_type;
32        typedef T * pointer;
33
34    public:
35        // TODO: think about a const pointer
36        typedef const pointer array_type;
37
38        typedef std::size_t size_type;
39        typedef std::ptrdiff_t difference_type;
40
41        typedef T value_type;
42        typedef const T  &const_reference;
43        typedef const T  *const_pointer;
44
45        typedef const_pointer const_iterator;
46        typedef std::reverse_iterator<const_iterator> const_reverse_iterator;
47
48        //
49        // typedefs required by vector concept
50        //
51
52        typedef dense_tag  storage_category;
53        typedef const vector_reference<const self_type>    const_closure_type;
54
55        c_array_view(size_type size, array_type data) :
56            size_(size), data_(data)
57        {}
58
59        ~c_array_view()
60        {}
61
62        //
63        // immutable methods of container concept
64        //
65
66        BOOST_UBLAS_INLINE
67        size_type size () const {
68            return size_;
69        }
70
71        BOOST_UBLAS_INLINE
72        const_reference operator [] (size_type i) const {
73            BOOST_UBLAS_CHECK (i < size_, bad_index ());
74            return data_ [i];
75        }
76
77        BOOST_UBLAS_INLINE
78        const_iterator begin () const {
79            return data_;
80        }
81        BOOST_UBLAS_INLINE
82        const_iterator end () const {
83            return data_ + size_;
84        }
85
86        BOOST_UBLAS_INLINE
87        const_reverse_iterator rbegin () const {
88            return const_reverse_iterator (end ());
89        }
90        BOOST_UBLAS_INLINE
91        const_reverse_iterator rend () const {
92            return const_reverse_iterator (begin ());
93        }
94
95    private:
96        size_type  size_;
97        array_type data_;
98    };
99
100
101    /** \brief Present existing arrays as compressed array based
102     *  sparse matrix.
103     *  This class provides CRS / CCS storage layout.
104     *
105     *  see also http://www.netlib.org/utk/papers/templates/node90.html
106     *
107     *       \param L layout type, either row_major or column_major
108     *       \param IB index base, use 0 for C indexing and 1 for
109     *       FORTRAN indexing of the internal index arrays. This
110     *       does not affect the operator()(int,int) where the first
111     *       row/column has always index 0.
112     *       \param IA index array type, e.g., int[]
113     *       \param TA value array type, e.g., double[]
114     */
115    template<class L, std::size_t IB, class IA, class JA, class TA>
116    class compressed_matrix_view:
117        public matrix_expression<compressed_matrix_view<L, IB, IA, JA, TA> > {
118
119    public:
120        typedef typename vector_view_traits<TA>::value_type value_type;
121
122    private:
123        typedef value_type &true_reference;
124        typedef value_type *pointer;
125        typedef const value_type *const_pointer;
126        typedef L layout_type;
127        typedef compressed_matrix_view<L, IB, IA, JA, TA> self_type;
128
129    public:
130#ifdef BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS
131        using matrix_expression<self_type>::operator ();
132#endif
133        // ISSUE require type consistency check
134        // is_convertable (IA::size_type, TA::size_type)
135        typedef typename boost::remove_cv<typename vector_view_traits<JA>::value_type>::type index_type;
136        // for compatibility, should be removed some day ...
137        typedef index_type size_type;
138        // size_type for the data arrays.
139        typedef typename vector_view_traits<JA>::size_type array_size_type;
140        typedef typename vector_view_traits<JA>::difference_type difference_type;
141        typedef const value_type & const_reference;
142
143        // do NOT define reference type, because class is read only
144        // typedef value_type & reference;
145
146        typedef IA rowptr_array_type;
147        typedef JA index_array_type;
148        typedef TA value_array_type;
149        typedef const matrix_reference<const self_type> const_closure_type;
150        typedef matrix_reference<self_type> closure_type;
151
152        // FIXME: define a corresponding temporary type
153        // typedef compressed_vector<T, IB, IA, TA> vector_temporary_type;
154
155        // FIXME: define a corresponding temporary type
156        // typedef self_type matrix_temporary_type;
157
158        typedef sparse_tag storage_category;
159        typedef typename L::orientation_category orientation_category;
160
161        //
162        // private types for internal use
163        //
164
165    private:
166        typedef typename vector_view_traits<index_array_type>::const_iterator const_subiterator_type;
167
168        //
169        // Construction and destruction
170        //
171    private:
172        /// private default constructor because data must be filled by caller
173        BOOST_UBLAS_INLINE
174        compressed_matrix_view () { }
175
176    public:
177        BOOST_UBLAS_INLINE
178        compressed_matrix_view (index_type n_rows, index_type n_cols, array_size_type nnz
179                                , const rowptr_array_type & iptr
180                                , const index_array_type & jptr
181                                , const value_array_type & values):
182            matrix_expression<self_type> (),
183            size1_ (n_rows), size2_ (n_cols), 
184            nnz_ (nnz),
185            index1_data_ (iptr), 
186            index2_data_ (jptr), 
187            value_data_ (values) {
188            storage_invariants ();
189        }
190
191        BOOST_UBLAS_INLINE
192        compressed_matrix_view(const compressed_matrix_view& o) :
193            size1_(size1_), size2_(size2_),
194            nnz_(nnz_),
195            index1_data_(index1_data_),
196            index2_data_(index2_data_),
197            value_data_(value_data_)
198        {}
199
200        //
201        // implement immutable iterator types
202        //
203
204        class const_iterator1 {};
205        class const_iterator2 {};
206
207        typedef reverse_iterator_base1<const_iterator1> const_reverse_iterator1;
208        typedef reverse_iterator_base2<const_iterator2> const_reverse_iterator2;
209
210        //
211        // implement all read only methods for the matrix expression concept
212        //
213
214        //! return the number of rows
215        index_type size1() const {
216            return size1_;
217        }
218
219        //! return the number of columns
220        index_type size2() const {
221            return size2_;
222        }
223
224        //! return value at position (i,j)
225        value_type operator()(index_type i, index_type j) const {
226            const_pointer p = find_element(i,j);
227            if (!p) {
228                return zero_;
229            } else {
230                return *p;
231            }
232        }
233       
234
235    private:
236        //
237        // private helper functions
238        //
239
240        const_pointer find_element (index_type i, index_type j) const {
241            index_type element1 (layout_type::index_M (i, j));
242            index_type element2 (layout_type::index_m (i, j));
243
244            const array_size_type itv      = zero_based( index1_data_[element1] );
245            const array_size_type itv_next = zero_based( index1_data_[element1+1] );
246
247            const_subiterator_type it_start = boost::next(vector_view_traits<index_array_type>::begin(index2_data_),itv);
248            const_subiterator_type it_end = boost::next(vector_view_traits<index_array_type>::begin(index2_data_),itv_next);
249            const_subiterator_type it = find_index_in_row(it_start, it_end, element2) ;
250           
251            if (it == it_end || *it != k_based (element2))
252                return 0;
253            return &value_data_ [it - vector_view_traits<index_array_type>::begin(index2_data_)];
254        }
255
256        const_subiterator_type find_index_in_row(const_subiterator_type it_start
257                                                 , const_subiterator_type it_end
258                                                 , index_type index) const {
259            return std::lower_bound( it_start
260                                     , it_end
261                                     , k_based (index) );
262        }
263
264
265    private:
266        void storage_invariants () const {
267            BOOST_UBLAS_CHECK (index1_data_ [layout_type::size_M (size1_, size2_)] == k_based (nnz_), external_logic ());
268        }
269       
270        index_type size1_;
271        index_type size2_;
272
273        array_size_type nnz_;
274
275        const rowptr_array_type & index1_data_;
276        const index_array_type & index2_data_;
277        const value_array_type & value_data_;
278
279        static const value_type zero_;
280
281        BOOST_UBLAS_INLINE
282        static index_type zero_based (index_type k_based_index) {
283            return k_based_index - IB;
284        }
285        BOOST_UBLAS_INLINE
286        static index_type k_based (index_type zero_based_index) {
287            return zero_based_index + IB;
288        }
289
290        friend class iterator1;
291        friend class iterator2;
292        friend class const_iterator1;
293        friend class const_iterator2;
294    };
295
296    template<class L, std::size_t IB, class IA, class JA, class TA  >
297    const typename compressed_matrix_view<L,IB,IA,JA,TA>::value_type
298    compressed_matrix_view<L,IB,IA,JA,TA>::zero_ = value_type/*zero*/();
299
300
301    template<class L, std::size_t IB, class IA, class JA, class TA  >
302    compressed_matrix_view<L,IB,IA,JA,TA>
303    make_compressed_matrix_view(typename vector_view_traits<JA>::value_type n_rows
304                                , typename vector_view_traits<JA>::value_type n_cols
305                                , typename vector_view_traits<JA>::size_type nnz
306                                , const IA & ia
307                                , const JA & ja
308                                , const TA & ta) {
309
310        return compressed_matrix_view<L,IB,IA,JA,TA>(n_rows, n_cols, nnz, ia, ja, ta);
311
312    }
313
314}}}
315
316#endif
Note: See TracBrowser for help on using the repository browser.