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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5588 - (show annotations)
Tue Apr 21 05:07:16 2015 UTC (5 years, 2 months ago) by caltinay
File MIME type: text/plain
File size: 11960 byte(s)
Update to copyright headers in cusp library.

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 /*
18 * Modifications to this file:
19 * Copyright (c) 2014-2015, The University of Queensland
20 * Licensed under the Apache License, Version 2.0.
21 *
22 */
23
24 /*! \file dia_matrix.h
25 * \brief Diagonal matrix format.
26 */
27
28 #pragma once
29
30 #include <cusp/detail/config.h>
31
32 #include <cusp/array1d.h>
33 #include <cusp/format.h>
34 #include <cusp/detail/matrix_base.h>
35 #include <cusp/detail/utils.h>
36
37 namespace cusp
38 {
39
40 /*! \addtogroup sparse_matrices Sparse Matrices
41 */
42
43 /*! \addtogroup sparse_matrix_containers Sparse Matrix Containers
44 * \ingroup sparse_matrices
45 * \{
46 */
47
48 // Forward definitions
49 struct column_major;
50 template<typename ValueType, class MemorySpace, class Orientation> class array2d;
51 template<typename Array, class Orientation> class array2d_view;
52 template <typename Array1, typename Array2, typename IndexType, typename ValueType, typename MemorySpace> class dia_matrix_view;
53
54 /*! \p dia_matrix : Diagonal matrix container
55 *
56 * \tparam IndexType Type used for matrix indices (e.g. \c int).
57 * \tparam ValueType Type used for matrix values (e.g. \c float).
58 * \tparam MemorySpace A memory space (e.g. \c cusp::host_memory or cusp::device_memory)
59 *
60 * \note The diagonal offsets should not contain duplicate entries.
61 *
62 * The following code snippet demonstrates how to create a 4-by-3
63 * \p dia_matrix on the host with 3 diagonals (6 total nonzeros)
64 * and then copies the matrix to the device.
65 *
66 * \code
67 * #include <cusp/dia_matrix.h>
68 * ...
69 *
70 * // allocate storage for (4,3) matrix with 6 nonzeros in 3 diagonals
71 * cusp::dia_matrix<int,float,cusp::host_memory> A(4,3,6,3);
72 *
73 * // initialize diagonal offsets
74 * A.diagonal_offsets[0] = -2;
75 * A.diagonal_offsets[1] = 0;
76 * A.diagonal_offsets[2] = 1;
77 *
78 * // initialize diagonal values
79 *
80 * // first diagonal
81 * A.values(0,2) = 0; // outside matrix
82 * A.values(1,2) = 0; // outside matrix
83 * A.values(2,0) = 40;
84 * A.values(3,0) = 60;
85 *
86 * // second diagonal
87 * A.values(0,1) = 10;
88 * A.values(1,1) = 0;
89 * A.values(2,1) = 50;
90 * A.values(3,1) = 50; // outside matrix
91 *
92 * // third diagonal
93 * A.values(0,2) = 20;
94 * A.values(1,2) = 30;
95 * A.values(2,2) = 0; // outside matrix
96 * A.values(3,2) = 0; // outside matrix
97 *
98 * // A now represents the following matrix
99 * // [10 20 0]
100 * // [ 0 0 30]
101 * // [40 0 50]
102 * // [ 0 60 0]
103 *
104 * // copy to the device
105 * cusp::dia_matrix<int,float,cusp::device_memory> B = A;
106 * \endcode
107 *
108 */
109 template <typename IndexType, typename ValueType, class MemorySpace>
110 class dia_matrix : public detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::dia_format>
111 {
112 typedef cusp::detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::dia_format> Parent;
113 public:
114 // TODO statically assert is_signed<IndexType>
115
116 /*! rebind matrix to a different MemorySpace
117 */
118 template<typename MemorySpace2>
119 struct rebind { typedef cusp::dia_matrix<IndexType, ValueType, MemorySpace2> type; };
120
121 /*! type of diagonal offsets array
122 */
123 typedef typename cusp::array1d<IndexType, MemorySpace> diagonal_offsets_array_type;
124
125 /*! type of values array
126 */
127 typedef typename cusp::array2d<ValueType, MemorySpace, cusp::column_major> values_array_type;
128
129 /*! equivalent container type
130 */
131 typedef typename cusp::dia_matrix<IndexType, ValueType, MemorySpace> container;
132
133 /*! equivalent view type
134 */
135 typedef typename cusp::dia_matrix_view<typename diagonal_offsets_array_type::view,
136 typename values_array_type::view,
137 IndexType, ValueType, MemorySpace> view;
138
139 /*! equivalent const_view type
140 */
141 typedef typename cusp::dia_matrix_view<typename diagonal_offsets_array_type::const_view,
142 typename values_array_type::const_view,
143 IndexType, ValueType, MemorySpace> const_view;
144
145 /*! Storage for the diagonal offsets.
146 */
147 diagonal_offsets_array_type diagonal_offsets;
148
149 /*! Storage for the nonzero entries of the DIA data structure.
150 */
151 values_array_type values;
152
153 /*! Indicator for matrix symmetry. If true, only non-negative diagonal
154 * offsets are stored.
155 */
156 bool symmetric;
157
158 /*! Construct an empty \p dia_matrix.
159 */
160 dia_matrix() : symmetric(false) {}
161
162 /*! Construct a \p dia_matrix with a specific shape, number of nonzero entries,
163 * and number of occupied diagonals.
164 *
165 * \param num_rows Number of rows.
166 * \param num_cols Number of columns.
167 * \param num_entries Number of nonzero matrix entries.
168 * \param num_diagonals Number of occupied diagonals.
169 * \param alignment Amount of padding used to align the data structure (default 32).
170 */
171 dia_matrix(size_t num_rows, size_t num_cols, size_t num_entries,
172 size_t num_diagonals, bool is_symmetric, size_t alignment = 32)
173 : Parent(num_rows, num_cols, num_entries),
174 diagonal_offsets(num_diagonals),
175 symmetric(is_symmetric)
176 {
177 // TODO use array2d constructor when it can accept pitch
178 values.resize(num_rows, num_diagonals, detail::round_up(num_rows, alignment));
179 }
180
181 /*! Construct a \p dia_matrix from another matrix.
182 *
183 * \param matrix Another sparse or dense matrix.
184 */
185 template <typename MatrixType>
186 dia_matrix(const MatrixType& matrix);
187
188 /*! Resize matrix dimensions and underlying storage
189 */
190 void resize(size_t num_rows, size_t num_cols, size_t num_entries,
191 size_t num_diagonals)
192 {
193 Parent::resize(num_rows, num_cols, num_entries);
194 diagonal_offsets.resize(num_diagonals);
195 values.resize(num_rows, num_diagonals);
196 }
197
198 /*! Resize matrix dimensions and underlying storage
199 */
200 void resize(size_t num_rows, size_t num_cols, size_t num_entries,
201 size_t num_diagonals, size_t alignment)
202 {
203 Parent::resize(num_rows, num_cols, num_entries);
204 diagonal_offsets.resize(num_diagonals);
205 values.resize(num_rows, num_diagonals, detail::round_up(num_rows, alignment));
206 }
207
208 /*! Swap the contents of two \p dia_matrix objects.
209 *
210 * \param matrix Another \p dia_matrix with the same IndexType and ValueType.
211 */
212 void swap(dia_matrix& matrix)
213 {
214 Parent::swap(matrix);
215 diagonal_offsets.swap(matrix.diagonal_offsets);
216 values.swap(matrix.values);
217 thrust::swap(symmetric, matrix.symmetric);
218 }
219
220 /*! Assignment from another matrix.
221 *
222 * \param matrix Another sparse or dense matrix.
223 */
224 template <typename MatrixType>
225 dia_matrix& operator=(const MatrixType& matrix);
226 }; // class dia_matrix
227 /*! \}
228 */
229
230 /*! \addtogroup sparse_matrix_views Sparse Matrix Views
231 * \ingroup sparse_matrices
232 * \{
233 */
234
235 /*! \p dia_matrix_view : Diagonal matrix view
236 *
237 * \tparam Array1 Type of \c diagonal_offsets
238 * \tparam Array2 Type of \c values array view
239 * \tparam IndexType Type used for matrix indices (e.g. \c int).
240 * \tparam ValueType Type used for matrix values (e.g. \c float).
241 * \tparam MemorySpace A memory space (e.g. \c cusp::host_memory or cusp::device_memory)
242 *
243 */
244 template <typename Array1,
245 typename Array2,
246 typename IndexType = typename Array1::value_type,
247 typename ValueType = typename Array2::value_type,
248 typename MemorySpace = typename cusp::minimum_space<typename Array1::memory_space, typename Array2::memory_space>::type >
249 class dia_matrix_view : public detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::dia_format>
250 {
251 typedef cusp::detail::matrix_base<IndexType,ValueType,MemorySpace,cusp::dia_format> Parent;
252 public:
253 /*! type of \c diagonal_offsets array
254 */
255 typedef Array1 diagonal_offsets_array_type;
256
257 /*! type of \c column_indices array
258 */
259 typedef Array2 values_array_type;
260
261 /*! equivalent container type
262 */
263 typedef typename cusp::dia_matrix<IndexType, ValueType, MemorySpace> container;
264
265 /*! equivalent view type
266 */
267 typedef typename cusp::dia_matrix_view<Array1, Array2, IndexType, ValueType, MemorySpace> view;
268
269 /*! Storage for the diagonal offsets.
270 */
271 diagonal_offsets_array_type diagonal_offsets;
272
273 /*! Storage for the nonzero entries of the DIA data structure.
274 */
275 values_array_type values;
276
277 /*! Construct an empty \p dia_matrix_view.
278 */
279 dia_matrix_view() {}
280
281 template <typename OtherArray1, typename OtherArray2>
282 dia_matrix_view(size_t num_rows, size_t num_cols, size_t num_entries,
283 OtherArray1& diagonal_offsets, OtherArray2& values)
284 : Parent(num_rows, num_cols, num_entries), diagonal_offsets(diagonal_offsets), values(values) {}
285
286 template <typename OtherArray1, typename OtherArray2>
287 dia_matrix_view(size_t num_rows, size_t num_cols, size_t num_entries,
288 const OtherArray1& diagonal_offsets, const OtherArray2& values)
289 : Parent(num_rows, num_cols, num_entries), diagonal_offsets(diagonal_offsets), values(values) {}
290
291 template <typename Matrix>
292 dia_matrix_view(Matrix& A)
293 : Parent(A), diagonal_offsets(A.diagonal_offsets), values(A.values) {}
294
295 template <typename Matrix>
296 dia_matrix_view(const Matrix& A)
297 : Parent(A), diagonal_offsets(A.diagonal_offsets), values(A.values) {}
298
299 /*! Resize matrix dimensions and underlying storage
300 */
301 void resize(size_t num_rows, size_t num_cols, size_t num_entries,
302 size_t num_diagonals)
303 {
304 Parent::resize(num_rows, num_cols, num_entries);
305 diagonal_offsets.resize(num_diagonals);
306 values.resize(num_rows, num_diagonals);
307 }
308
309 /*! Resize matrix dimensions and underlying storage
310 */
311 void resize(size_t num_rows, size_t num_cols, size_t num_entries,
312 size_t num_diagonals, size_t alignment)
313 {
314 Parent::resize(num_rows, num_cols, num_entries);
315 diagonal_offsets.resize(num_diagonals);
316 values.resize(num_rows, num_diagonals, detail::round_up(num_rows, alignment));
317 }
318 }; // class dia_matrix_view
319
320
321 template <typename Array1,
322 typename Array2>
323 dia_matrix_view<Array1,Array2>
324 make_dia_matrix_view(size_t num_rows,
325 size_t num_cols,
326 size_t num_entries,
327 Array1 diagonal_offsets,
328 Array2 values);
329
330 template <typename Array1,
331 typename Array2,
332 typename IndexType,
333 typename ValueType,
334 typename MemorySpace>
335 dia_matrix_view<Array1,Array2,IndexType,ValueType,MemorySpace>
336 make_dia_matrix_view(const dia_matrix_view<Array1,Array2,IndexType,ValueType,MemorySpace>& m);
337
338 template <typename IndexType, typename ValueType, class MemorySpace>
339 typename dia_matrix<IndexType,ValueType,MemorySpace>::view
340 make_dia_matrix_view(dia_matrix<IndexType,ValueType,MemorySpace>& m);
341
342 template <typename IndexType, typename ValueType, class MemorySpace>
343 typename dia_matrix<IndexType,ValueType,MemorySpace>::const_view
344 make_dia_matrix_view(const dia_matrix<IndexType,ValueType,MemorySpace>& m);
345 /*! \} // end Views
346 */
347
348 } // end namespace cusp
349
350 #include <cusp/array2d.h>
351 #include <cusp/detail/dia_matrix.inl>
352

  ViewVC Help
Powered by ViewVC 1.1.26