/[escript]/trunk/cusplibrary/cusp/hyb_matrix.h
ViewVC logotype

Contents of /trunk/cusplibrary/cusp/hyb_matrix.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5148 - (show annotations)
Mon Sep 15 01:25:23 2014 UTC (5 years, 8 months ago) by caltinay
File MIME type: text/plain
File size: 11685 byte(s)
Merging ripley diagonal storage + CUDA support into trunk.
Options file version has been incremented due to new options
'cuda' and 'nvccflags'.

1 /*
2 * Copyright 2008-2009 NVIDIA Corporation
3 *
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
7 *
8 * http://www.apache.org/licenses/LICENSE-2.0
9 *
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13 * See the License for the specific language governing permissions and
14 * limitations under the License.
15 */
16
17 /*! \file hyb_matrix.h
18 * \brief Hybrid ELL/COO matrix format
19 */
20
21 #pragma once
22
23 #include <cusp/detail/config.h>
24
25 #include <cusp/array1d.h>
26 #include <cusp/format.h>
27 #include <cusp/detail/matrix_base.h>
28
29 namespace cusp
30 {
31 // Forward definitions
32 template <typename IndexType, typename ValueType, class MemorySpace> class ell_matrix;
33 template <typename IndexType, typename ValueType, class MemorySpace> class coo_matrix;
34 template <typename Matrix1, typename Matrix2, typename IndexType, typename ValueType, class MemorySpace> class hyb_matrix_view;
35
36 /*! \addtogroup sparse_matrices Sparse Matrices
37 */
38
39 /*! \addtogroup sparse_matrix_containers Sparse Matrix Containers
40 * \ingroup sparse_matrices
41 * \{
42 */
43
44 /*! \p hyb_matrix : Hybrid ELL/COO matrix container
45 *
46 * The \p hyb_matrix is a combination of the \p ell_matrix and
47 * \p coo_matrix formats. Specifically, the \p hyb_matrix format
48 * splits a matrix into two portions, one stored in ELL format
49 * and one stored in COO format.
50 *
51 * While the ELL format is well-suited to vector and SIMD
52 * architectures, its efficiency rapidly degrades when the number of
53 * nonzeros per matrix row varies. In contrast, the storage efficiency of
54 * the COO format is invariant to the distribution of nonzeros per row, and
55 * the use of segmented reduction makes its performance largely invariant
56 * as well. To obtain the advantages of both, we combine these
57 * into a hybrid ELL/COO format.
58 *
59 * The purpose of the HYB format is to store the typical number of
60 * nonzeros per row in the ELL data structure and the remaining entries of
61 * exceptional rows in the COO format.
62 *
63 * \tparam IndexType Type used for matrix indices (e.g. \c int).
64 * \tparam ValueType Type used for matrix values (e.g. \c float).
65 * \tparam MemorySpace A memory space (e.g. \c cusp::host_memory or cusp::device_memory)
66 *
67 * \note The \p ell_matrix entries must be sorted by column index.
68 * \note The \p ell_matrix entries within each row should be shifted to the left.
69 * \note The \p coo_matrix entries must be sorted by row index.
70 * \note The matrix should not contain duplicate entries.
71 *
72 * The following code snippet demonstrates how to create a \p hyb_matrix.
73 * In practice we usually do not construct the HYB format directly and
74 * instead convert from a simpler format such as (COO, CSR) into HYB.
75 *
76 * \code
77 * #include <cusp/hyb_matrix.h>
78 * ...
79 *
80 * // allocate storage for (4,3) matrix with 8 nonzeros
81 * // ELL portion has 5 nonzeros and storage for 2 nonzeros per row
82 * // COO portion has 3 nonzeros
83 *
84 * cusp::hyb_matrix<int, float, cusp::host_memory> A(3, 4, 5, 3, 2);
85 *
86 * // Initialize A to represent the following matrix
87 * // [10 20 30 40]
88 * // [ 0 50 0 0]
89 * // [60 0 70 80]
90 *
91 * // A is split into ELL and COO parts as follows
92 * // [10 20 30 40] [10 20 0 0] [ 0 0 30 40]
93 * // [ 0 50 0 0] = [ 0 50 0 0] + [ 0 0 0 0]
94 * // [60 0 70 80] [60 0 70 0] [ 0 0 0 80]
95 *
96 *
97 * // Initialize ELL part
98 *
99 * // X is used to fill unused entries in the ELL portion of the matrix
100 * const int X = cusp::ell_matrix<int,float,cusp::host_memory>::invalid_index;
101 *
102 * // first row
103 * A.ell.column_indices(0,0) = 0; A.ell.values(0,0) = 10;
104 * A.ell.column_indices(0,1) = 1; A.ell.values(0,1) = 20;
105 *
106 * // second row
107 * A.ell.column_indices(1,0) = 1; A.ell.values(1,0) = 50; // shifted to leftmost position
108 * A.ell.column_indices(1,1) = X; A.ell.values(1,1) = 0; // padding
109 *
110 * // third row
111 * A.ell.column_indices(2,0) = 0; A.ell.values(2,0) = 60;
112 * A.ell.column_indices(2,1) = 2; A.ell.values(2,1) = 70; // shifted to leftmost position
113 *
114 *
115 * // Initialize COO part
116 * A.coo.row_indices[0] = 0; A.coo.column_indices[0] = 2; A.coo.values[0] = 30;
117 * A.coo.row_indices[1] = 0; A.coo.column_indices[1] = 3; A.coo.values[1] = 40;
118 * A.coo.row_indices[2] = 2; A.coo.column_indices[2] = 3; A.coo.values[2] = 80;
119 *
120 * \endcode
121 *
122 * \see \p ell_matrix
123 * \see \p coo_matrix
124 */
125 template <typename IndexType, typename ValueType, class MemorySpace>
126 class hyb_matrix : public detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::hyb_format>
127 {
128 typedef cusp::detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::hyb_format> Parent;
129 public:
130 /*! rebind matrix to a different MemorySpace
131 */
132 template<typename MemorySpace2>
133 struct rebind { typedef cusp::hyb_matrix<IndexType, ValueType, MemorySpace2> type; };
134
135 /*! equivalent container type
136 */
137 typedef typename cusp::hyb_matrix<IndexType, ValueType, MemorySpace> container;
138
139 /*! equivalent view type
140 */
141 typedef typename cusp::hyb_matrix_view<typename cusp::ell_matrix<IndexType,ValueType,MemorySpace>::view,
142 typename cusp::coo_matrix<IndexType,ValueType,MemorySpace>::view,
143 IndexType, ValueType, MemorySpace> view;
144
145 /*! equivalent const_view type
146 */
147 typedef typename cusp::hyb_matrix_view<typename cusp::ell_matrix<IndexType,ValueType,MemorySpace>::const_view,
148 typename cusp::coo_matrix<IndexType,ValueType,MemorySpace>::const_view,
149 IndexType, ValueType, MemorySpace> const_view;
150
151 /*! type of \p ELL portion of the HYB structure
152 */
153 typedef cusp::ell_matrix<IndexType,ValueType,MemorySpace> ell_matrix_type;
154
155 /*! type of \p COO portion of the HYB structure
156 */
157 typedef cusp::coo_matrix<IndexType,ValueType,MemorySpace> coo_matrix_type;
158
159 /*! Storage for the \p ell_matrix portion.
160 */
161 ell_matrix_type ell;
162
163 /*! Storage for the \p ell_matrix portion.
164 */
165 coo_matrix_type coo;
166
167 /*! Construct an empty \p hyb_matrix.
168 */
169 hyb_matrix() {}
170
171 /*! Construct a \p hyb_matrix with a specific shape and separation into ELL and COO portions.
172 *
173 * \param num_rows Number of rows.
174 * \param num_cols Number of columns.
175 * \param num_ell_entries Number of nonzero matrix entries in the ELL portion.
176 * \param num_coo_entries Number of nonzero matrix entries in the ELL portion.
177 * \param num_entries_per_row Maximum number of nonzeros per row in the ELL portion.
178 * \param alignment Amount of padding used to align the ELL data structure (default 32).
179 */
180 hyb_matrix(IndexType num_rows, IndexType num_cols,
181 IndexType num_ell_entries, IndexType num_coo_entries,
182 IndexType num_entries_per_row, IndexType alignment = 32)
183 : Parent(num_rows, num_cols, num_ell_entries + num_coo_entries),
184 ell(num_rows, num_cols, num_ell_entries, num_entries_per_row, alignment),
185 coo(num_rows, num_cols, num_coo_entries) {}
186
187 // TODO remove default alignment of 32
188
189 /*! Construct a \p hyb_matrix from another matrix.
190 *
191 * \param matrix Another sparse or dense matrix.
192 */
193 template <typename MatrixType>
194 hyb_matrix(const MatrixType& matrix);
195
196 /*! Resize matrix dimensions and underlying storage
197 */
198 void resize(IndexType num_rows, IndexType num_cols,
199 IndexType num_ell_entries, IndexType num_coo_entries,
200 IndexType num_entries_per_row, IndexType alignment = 32)
201 {
202 Parent::resize(num_rows, num_cols, num_ell_entries + num_coo_entries);
203 ell.resize(num_rows, num_cols, num_ell_entries, num_entries_per_row, alignment);
204 coo.resize(num_rows, num_cols, num_coo_entries);
205 }
206
207 /*! Swap the contents of two \p hyb_matrix objects.
208 *
209 * \param matrix Another \p hyb_matrix with the same IndexType and ValueType.
210 */
211 void swap(hyb_matrix& matrix)
212 {
213 Parent::swap(matrix);
214 ell.swap(matrix.ell);
215 coo.swap(matrix.coo);
216 }
217
218 /*! Assignment from another matrix.
219 *
220 * \param matrix Another sparse or dense matrix.
221 */
222 template <typename MatrixType>
223 hyb_matrix& operator=(const MatrixType& matrix);
224 }; // class hyb_matrix
225 /*! \}
226 */
227
228
229 /*! \addtogroup sparse_matrix_views Sparse Matrix Views
230 * \ingroup sparse_matrices
231 * \{
232 */
233
234 /*! \p hyb_matrix_view : Hybrid ELL/COO matrix view
235 *
236 * \tparam Matrix1 Type of \c ell
237 * \tparam Matrix2 Type of \c coo
238 * \tparam IndexType Type used for matrix indices (e.g. \c int).
239 * \tparam ValueType Type used for matrix values (e.g. \c float).
240 * \tparam MemorySpace A memory space (e.g. \c cusp::host_memory or cusp::device_memory)
241 *
242 */
243 template <typename Matrix1,
244 typename Matrix2,
245 typename IndexType = typename Matrix1::index_type,
246 typename ValueType = typename Matrix1::value_type,
247 typename MemorySpace = typename cusp::minimum_space<typename Matrix1::memory_space, typename Matrix2::memory_space>::type >
248 class hyb_matrix_view : public detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::hyb_format>
249 {
250 typedef cusp::detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::hyb_format> Parent;
251 public:
252 /*! type of \p ELL portion of the HYB structure
253 */
254 typedef Matrix1 ell_matrix_type;
255
256 /*! type of \p COO portion of the HYB structure
257 */
258 typedef Matrix2 coo_matrix_type;
259
260 /*! equivalent container type
261 */
262 typedef typename cusp::hyb_matrix<IndexType, ValueType, MemorySpace> container;
263
264 /*! equivalent view type
265 */
266 typedef typename cusp::hyb_matrix_view<Matrix1, Matrix2, IndexType, ValueType, MemorySpace> view;
267
268 /*! View to the \p ELL portion of the HYB structure.
269 */
270 ell_matrix_type ell;
271
272 /*! View to the \p COO portion of the HYB structure.
273 */
274 coo_matrix_type coo;
275
276 /*! Construct an empty \p hyb_matrix_view.
277 */
278 hyb_matrix_view() {}
279
280 template <typename OtherMatrix1, typename OtherMatrix2>
281 hyb_matrix_view(OtherMatrix1& ell, OtherMatrix2& coo)
282 : Parent(ell.num_rows, ell.num_cols, ell.num_entries + coo.num_entries), ell(ell), coo(coo) {}
283
284 template <typename OtherMatrix1, typename OtherMatrix2>
285 hyb_matrix_view(const OtherMatrix1& ell, const OtherMatrix2& coo)
286 : Parent(ell.num_rows, ell.num_cols, ell.num_entries + coo.num_entries), ell(ell), coo(coo) {}
287
288 template <typename Matrix>
289 hyb_matrix_view(Matrix& A)
290 : Parent(A), ell(A.ell), coo(A.coo) {}
291
292 template <typename Matrix>
293 hyb_matrix_view(const Matrix& A)
294 : Parent(A), ell(A.ell), coo(A.coo) {}
295
296 /*! Resize matrix dimensions and underlying storage
297 */
298 void resize(size_t num_rows, size_t num_cols,
299 size_t num_ell_entries, size_t num_coo_entries,
300 size_t num_entries_per_row, size_t alignment = 32)
301 {
302 Parent::resize(num_rows, num_cols, num_ell_entries + num_coo_entries);
303 ell.resize(num_rows, num_cols, num_ell_entries, num_entries_per_row, alignment);
304 coo.resize(num_rows, num_cols, num_coo_entries);
305 }
306 };
307 /*! \} // end Views
308 */
309
310 } // end namespace cusp
311
312 #include <cusp/detail/hyb_matrix.inl>
313

  ViewVC Help
Powered by ViewVC 1.1.26