/[escript]/branches/diaplayground/cusplibrary/cusp/array2d.h
ViewVC logotype

Contents of /branches/diaplayground/cusplibrary/cusp/array2d.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4956 - (show annotations)
Tue May 20 12:34:41 2014 UTC (6 years ago) by caltinay
File MIME type: text/plain
File size: 17960 byte(s)
minor changes to cusp to allow compilation with -pedantic and -Werror.
One bug fix in spmv_dia.

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 array2d.h
18 * \brief Two-dimensional array
19 */
20
21 #pragma once
22
23 #include <cusp/detail/config.h>
24
25 #include <cusp/memory.h>
26 #include <cusp/format.h>
27 #include <cusp/array1d.h>
28 #include <cusp/detail/matrix_base.h>
29 #include <cusp/detail/utils.h>
30
31 #include <thrust/functional.h>
32
33 // TODO move indexing stuff to /detail/
34
35 namespace cusp
36 {
37
38 struct row_major {};
39 struct column_major {};
40
41 // forward definitions
42 template<typename Array, class Orientation> class array2d_view;
43
44 namespace detail
45 {
46 // (i,j) -> major dimension
47 // (i,j) -> minor dimension
48 // logical n -> (i,j)
49 // (i,j) -> logical n
50 // (i,j) -> physical n
51 // logical n -> physical n
52 // logical n -> physical n (translated)
53
54 template <typename IndexType>
55 __host__ __device__
56 IndexType minor_dimension(IndexType num_rows, IndexType num_cols, row_major) {
57 return num_cols;
58 }
59
60 template <typename IndexType>
61 __host__ __device__
62 IndexType minor_dimension(IndexType num_rows, IndexType num_cols, column_major) {
63 return num_rows;
64 }
65
66 template <typename IndexType>
67 __host__ __device__
68 IndexType major_dimension(IndexType num_rows, IndexType num_cols, row_major) {
69 return num_rows;
70 }
71
72 template <typename IndexType>
73 __host__ __device__
74 IndexType major_dimension(IndexType num_rows, IndexType num_cols, column_major) {
75 return num_cols;
76 }
77
78 // convert logical linear index into a logical (i,j) index
79 template <typename IndexType>
80 __host__ __device__
81 IndexType linear_index_to_row_index(IndexType linear_index, IndexType num_rows, IndexType num_cols, row_major) {
82 return linear_index / num_cols;
83 }
84
85 template <typename IndexType>
86 __host__ __device__
87 IndexType linear_index_to_col_index(IndexType linear_index, IndexType num_rows, IndexType num_cols, row_major) {
88 return linear_index % num_cols;
89 }
90
91 template <typename IndexType>
92 __host__ __device__
93 IndexType linear_index_to_row_index(IndexType linear_index, IndexType num_rows, IndexType num_cols, column_major) {
94 return linear_index % num_rows;
95 }
96
97 template <typename IndexType>
98 __host__ __device__
99 IndexType linear_index_to_col_index(IndexType linear_index, IndexType num_rows, IndexType num_cols, column_major) {
100 return linear_index / num_rows;
101 }
102
103 // convert a logical (i,j) index into a physical linear index
104 template <typename IndexType>
105 __host__ __device__
106 IndexType index_of(IndexType i, IndexType j, IndexType pitch, row_major) {
107 return i * pitch + j;
108 }
109
110 template <typename IndexType>
111 __host__ __device__
112 IndexType index_of(IndexType i, IndexType j, IndexType pitch, column_major) {
113 return j * pitch + i;
114 }
115
116 template <typename IndexType, typename Orientation>
117 __host__ __device__
118 IndexType logical_to_physical(IndexType linear_index, IndexType num_rows, IndexType num_cols, IndexType pitch, Orientation)
119 {
120 IndexType i = linear_index_to_row_index(linear_index, num_rows, num_cols, Orientation());
121 IndexType j = linear_index_to_col_index(linear_index, num_rows, num_cols, Orientation());
122
123 return index_of(i, j, pitch, Orientation());
124 }
125
126 // convert logical linear index in the source into a physical linear index in the destination
127 template <typename IndexType, typename Orientation1, typename Orientation2>
128 __host__ __device__
129 IndexType logical_to_other_physical(IndexType linear_index, IndexType num_rows, IndexType num_cols, IndexType pitch, Orientation1, Orientation2)
130 {
131 IndexType i = linear_index_to_row_index(linear_index, num_rows, num_cols, Orientation1());
132 IndexType j = linear_index_to_col_index(linear_index, num_rows, num_cols, Orientation1());
133
134 return index_of(i, j, pitch, Orientation2());
135 }
136
137 // functors
138 template <typename IndexType, typename Orientation>
139 struct logical_to_physical_functor : public thrust::unary_function<IndexType,IndexType>
140 {
141 IndexType num_rows, num_cols, pitch;
142
143 logical_to_physical_functor(IndexType num_rows, IndexType num_cols, IndexType pitch)
144 : num_rows(num_rows), num_cols(num_cols), pitch(pitch) {}
145
146 __host__ __device__
147 IndexType operator()(const IndexType i) const
148 {
149 return logical_to_physical(i, num_rows, num_cols, pitch, Orientation());
150 }
151 };
152
153 template <typename IndexType, typename Orientation1, typename Orientation2>
154 struct logical_to_other_physical_functor : public thrust::unary_function<IndexType,IndexType>
155 {
156 IndexType num_rows, num_cols, pitch;
157
158 logical_to_other_physical_functor(IndexType num_rows, IndexType num_cols, IndexType pitch)
159 : num_rows(num_rows), num_cols(num_cols), pitch(pitch) {}
160
161 __host__ __device__
162 IndexType operator()(const IndexType i) const
163 {
164 return logical_to_other_physical(i, num_rows, num_cols, pitch, Orientation1(), Orientation2());
165 }
166 };
167
168 template <typename Iterator, bool same_orientation>
169 struct row_or_column_view {};
170
171 template <typename Iterator>
172 struct row_or_column_view<Iterator,true>
173 {
174 typedef cusp::array1d_view<Iterator> ArrayType;
175
176 template< typename Array >
177 static ArrayType get_array(Array& A, size_t i ) {
178 return ArrayType(A.values.begin() + A.pitch * i,
179 A.values.begin() + A.pitch * i + cusp::detail::minor_dimension(A.num_rows, A.num_cols, typename Array::orientation()));
180 }
181 };
182
183 template <typename Iterator>
184 struct row_or_column_view<Iterator,false>
185 {
186 typedef typename cusp::detail::strided_range<Iterator> StrideType;
187 typedef cusp::array1d_view<typename StrideType::iterator> ArrayType;
188
189 template< typename Array >
190 static ArrayType get_array(Array& A, size_t i ) {
191 cusp::detail::strided_range<Iterator> strided_range(A.values.begin() + i, A.values.end(), A.pitch);
192 return ArrayType(strided_range.begin(), strided_range.end());
193 }
194 };
195
196 } // end namespace detail
197
198 // TODO document mapping of (i,j) onto values[pitch * i + j] or values[pitch * j + i]
199 // TODO document that array2d operations will try to respect .pitch of destination argument
200
201 /*! \addtogroup arrays Arrays
202 */
203
204 /*! \addtogroup array_containers Array Containers
205 * \ingroup arrays
206 * \{
207 */
208
209 /*! \p array2d : One-dimensional array container
210 *
211 * \tparam T value_type of the array
212 * \tparam MemorySpace memory space of the array (cusp::host_memory or cusp::device_memory)
213 * \tparam Orientation orientation of the array (cusp::row_major or cusp::column_major)
214 *
215 * \TODO example
216 */
217 template<typename ValueType, class MemorySpace, class Orientation = cusp::row_major>
218 class array2d : public cusp::detail::matrix_base<int,ValueType,MemorySpace,cusp::array2d_format>
219 {
220 typedef typename cusp::detail::matrix_base<int,ValueType,MemorySpace,cusp::array2d_format> Parent;
221
222 public:
223 typedef Orientation orientation;
224
225 template<typename MemorySpace2>
226 struct rebind {
227 typedef cusp::array2d<ValueType, MemorySpace2, Orientation> type;
228 };
229
230 typedef typename cusp::array1d<ValueType, MemorySpace> values_array_type;
231
232 /*! equivalent container type
233 */
234 typedef typename cusp::array2d<ValueType, MemorySpace, Orientation> container;
235
236 /*! equivalent view type
237 */
238 typedef typename cusp::array2d_view<typename values_array_type::view, Orientation> view;
239
240 /*! equivalent const_view type
241 */
242 typedef typename cusp::array2d_view<typename values_array_type::const_view, Orientation> const_view;
243
244 /*! array1d_view of a single row
245 */
246 typedef cusp::detail::row_or_column_view<typename values_array_type::iterator,thrust::detail::is_same<Orientation,cusp::row_major>::value> row_view_type;
247 typedef typename row_view_type::ArrayType row_view;
248
249 /*! array1d_view of a single column
250 */
251 typedef cusp::detail::row_or_column_view<typename values_array_type::iterator,thrust::detail::is_same<Orientation,cusp::column_major>::value> column_view_type;
252 typedef typename column_view_type::ArrayType column_view;
253
254 /*! const array1d_view of a single row
255 */
256 typedef cusp::detail::row_or_column_view<typename values_array_type::const_iterator,thrust::detail::is_same<Orientation,cusp::row_major>::value> const_row_view_type;
257 typedef typename const_row_view_type::ArrayType const_row_view;
258
259 /*! const array1d_view of a single column
260 */
261 typedef cusp::detail::row_or_column_view<typename values_array_type::const_iterator,thrust::detail::is_same<Orientation,cusp::column_major>::value> const_column_view_type;
262 typedef typename const_column_view_type::ArrayType const_column_view;
263
264 // minor_dimension + padding
265 size_t pitch;
266
267 values_array_type values;
268
269 // construct empty matrix
270 array2d()
271 : Parent(), pitch(0), values(0) {}
272
273 // construct matrix with given shape and number of entries
274 array2d(size_t num_rows, size_t num_cols)
275 : Parent(num_rows, num_cols, num_rows * num_cols),
276 pitch(cusp::detail::minor_dimension(num_rows, num_cols, orientation())),
277 values(num_rows * num_cols) {}
278
279 // construct matrix with given shape, number of entries and fill with a given value
280 array2d(size_t num_rows, size_t num_cols, const ValueType& value)
281 : Parent(num_rows, num_cols, num_rows * num_cols),
282 pitch(cusp::detail::minor_dimension(num_rows, num_cols, orientation())),
283 values(num_rows * num_cols, value) {}
284
285 // construct matrix with given shape, number of entries, pitch and fill with a given value
286 array2d(size_t num_rows, size_t num_cols, const ValueType& value, size_t pitch)
287 : Parent(num_rows, num_cols, num_rows * num_cols),
288 pitch(pitch),
289 values(pitch * cusp::detail::major_dimension(num_rows, num_cols, orientation()), value)
290 {
291 if (pitch < cusp::detail::minor_dimension(num_rows, num_cols, orientation()))
292 throw cusp::invalid_input_exception("pitch cannot be less than minor dimension");
293 }
294
295 // construct from another matrix
296 template <typename MatrixType>
297 array2d(const MatrixType& matrix);
298
299 typename values_array_type::reference operator()(const size_t i, const size_t j)
300 {
301 return values[cusp::detail::index_of(i, j, pitch, orientation())];
302 }
303
304 typename values_array_type::const_reference operator()(const size_t i, const size_t j) const
305 {
306 return values[cusp::detail::index_of(i, j, pitch, orientation())];
307 }
308
309 void resize(size_t num_rows, size_t num_cols, size_t pitch)
310 {
311 if (pitch < cusp::detail::minor_dimension(num_rows, num_cols, orientation()))
312 throw cusp::invalid_input_exception("pitch cannot be less than minor dimension");
313
314 values.resize(pitch * cusp::detail::major_dimension(num_rows, num_cols, orientation()));
315
316 this->num_rows = num_rows;
317 this->num_cols = num_cols;
318 this->pitch = pitch;
319 this->num_entries = num_rows * num_cols;
320 }
321
322 void resize(size_t num_rows, size_t num_cols)
323 {
324 // preserve .pitch if possible
325 if (this->num_rows == num_rows && this->num_cols == num_cols)
326 return;
327
328 resize(num_rows, num_cols, cusp::detail::minor_dimension(num_rows, num_cols, orientation()));
329 }
330
331 void swap(array2d& matrix)
332 {
333 Parent::swap(matrix);
334 thrust::swap(this->pitch, matrix.pitch);
335 values.swap(matrix.values);
336 }
337
338 row_view row(size_t i)
339 {
340 return row_view_type::get_array(*this, i);
341 }
342
343 column_view column(size_t i)
344 {
345 return column_view_type::get_array(*this, i);
346 }
347
348 const_row_view row(size_t i) const
349 {
350 return const_row_view_type::get_array(*this, i);
351 }
352
353 const_column_view column(size_t i) const
354 {
355 return const_column_view_type::get_array(*this, i);
356 }
357
358 array2d& operator=(const array2d& matrix);
359
360 template <typename MatrixType>
361 array2d& operator=(const MatrixType& matrix);
362
363 }; // class array2d
364 /*! \}
365 */
366
367 /*! \addtogroup array_views Array Views
368 * \ingroup arrays
369 * \{
370 */
371
372 /*! \p array2d_view : One-dimensional array view
373 *
374 * \tparam Array Underlying one-dimensional array view
375 * \tparam Orientation orientation of the array (cusp::row_major or cusp::column_major)
376 *
377 * \TODO example
378 */
379 template<typename Array, class Orientation = cusp::row_major>
380 class array2d_view : public cusp::detail::matrix_base<int, typename Array::value_type,typename Array::memory_space, cusp::array2d_format>
381 {
382 typedef cusp::detail::matrix_base<int, typename Array::value_type,typename Array::memory_space, cusp::array2d_format> Parent;
383 public:
384 typedef Orientation orientation;
385
386 typedef Array values_array_type;
387
388 values_array_type values;
389
390 /*! equivalent container type
391 */
392 typedef typename cusp::array2d<typename Parent::value_type, typename Parent::memory_space, Orientation> container;
393
394 /*! equivalent view type
395 */
396 typedef typename cusp::array2d_view<Array, Orientation> view;
397
398 /*! array1d_view of a single row
399 */
400 typedef cusp::detail::row_or_column_view<typename values_array_type::iterator,thrust::detail::is_same<Orientation,cusp::row_major>::value> row_view_type;
401 typedef typename row_view_type::ArrayType row_view;
402
403 /*! array1d_view of a single column
404 */
405 typedef cusp::detail::row_or_column_view<typename values_array_type::iterator,thrust::detail::is_same<Orientation,cusp::column_major>::value> column_view_type;
406 typedef typename column_view_type::ArrayType column_view;
407
408 // minor_dimension + padding
409 size_t pitch;
410
411 // construct empty view
412 array2d_view(void)
413 : Parent(), values(0), pitch(0) {}
414
415 array2d_view(const array2d_view& a)
416 : Parent(a), values(a.values), pitch(a.pitch) {}
417
418 // TODO handle different Orientation (pitch = major)
419 //template <typename Array2, typename Orientation2>
420 //array2d_view(const array2d_view<Array2,Orientation2>& A)
421
422 // TODO check values.size()
423
424 // construct from array2d container
425 array2d_view( array2d<typename Parent::value_type, typename Parent::memory_space, orientation>& a)
426 : Parent(a), values(a.values), pitch(a.pitch) {}
427
428 array2d_view(const array2d<typename Parent::value_type, typename Parent::memory_space, orientation>& a)
429 : Parent(a), values(a.values), pitch(a.pitch) {}
430
431 template <typename Array2>
432 array2d_view(size_t num_rows, size_t num_cols, size_t pitch, Array2& values)
433 : Parent(num_rows, num_cols, num_rows * num_cols), pitch(pitch), values(values) {}
434
435 template <typename Array2>
436 array2d_view(size_t num_rows, size_t num_cols, size_t pitch, const Array2& values)
437 : Parent(num_rows, num_cols, num_rows * num_cols), pitch(pitch), values(values) {}
438
439 typename values_array_type::reference operator()(const size_t i, const size_t j) const
440 {
441 return values[detail::index_of(i, j, pitch, orientation())];
442 }
443
444 void resize(size_t num_rows, size_t num_cols, size_t pitch)
445 {
446 if (pitch < cusp::detail::minor_dimension(num_rows, num_cols, orientation()))
447 throw cusp::invalid_input_exception("pitch cannot be less than minor dimension");
448
449 values.resize(pitch * cusp::detail::major_dimension(num_rows, num_cols, orientation()));
450
451 this->num_rows = num_rows;
452 this->num_cols = num_cols;
453 this->pitch = pitch;
454 this->num_entries = num_rows * num_cols;
455 }
456
457 void resize(size_t num_rows, size_t num_cols)
458 {
459 // preserve .pitch if possible
460 if (this->num_rows == num_rows && this->num_cols == num_cols)
461 return;
462
463 resize(num_rows, num_cols, cusp::detail::minor_dimension(num_rows, num_cols, orientation()));
464 }
465
466 row_view row(size_t i)
467 {
468 return row_view_type::get_array(*this, i);
469 }
470
471 column_view column(size_t i)
472 {
473 return column_view_type::get_array(*this, i);
474 }
475
476 row_view row(size_t i) const
477 {
478 return row_view_type::get_array(*this, i);
479 }
480
481 column_view column(size_t i) const
482 {
483 return column_view_type::get_array(*this, i);
484 }
485 }; // class array2d_view
486
487
488
489 template <typename Iterator, typename Orientation>
490 array2d_view<typename cusp::array1d_view<Iterator>,Orientation>
491 make_array2d_view(size_t num_rows, size_t num_cols, size_t pitch, const cusp::array1d_view<Iterator>& values, Orientation)
492 {
493 return array2d_view<typename cusp::array1d_view<Iterator>,Orientation>(num_rows, num_cols, pitch, values);
494 }
495
496 template <typename Array, typename Orientation>
497 array2d_view<Array,Orientation>
498 make_array2d_view(const array2d_view<Array, Orientation>& a)
499 {
500 return array2d_view<Array,Orientation>(a);
501 }
502
503 template<typename ValueType, class MemorySpace, class Orientation>
504 array2d_view<typename cusp::array1d_view<typename cusp::array1d<ValueType,MemorySpace>::iterator >, Orientation>
505 make_array2d_view(cusp::array2d<ValueType,MemorySpace,Orientation>& a)
506 {
507 return cusp::make_array2d_view(a.num_rows, a.num_cols, a.pitch, cusp::make_array1d_view(a.values), Orientation());
508 }
509
510 template<typename ValueType, class MemorySpace, class Orientation>
511 array2d_view<typename cusp::array1d_view<typename cusp::array1d<ValueType,MemorySpace>::const_iterator >, Orientation>
512 make_array2d_view(const cusp::array2d<ValueType,MemorySpace,Orientation>& a)
513 {
514 return cusp::make_array2d_view(a.num_rows, a.num_cols, a.pitch, cusp::make_array1d_view(a.values), Orientation());
515 }
516 /*! \}
517 */
518
519 } // end namespace cusp
520
521 #include <cusp/detail/array2d.inl>
522

  ViewVC Help
Powered by ViewVC 1.1.26