ViennaCL - The Vienna Computing Library  1.5.1
misc_operations.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_CUDA_MISC_OPERATIONS_HPP_
2 #define VIENNACL_LINALG_CUDA_MISC_OPERATIONS_HPP_
3 
4 /* =========================================================================
5  Copyright (c) 2010-2014, Institute for Microelectronics,
6  Institute for Analysis and Scientific Computing,
7  TU Wien.
8  Portions of this software are copyright by UChicago Argonne, LLC.
9 
10  -----------------
11  ViennaCL - The Vienna Computing Library
12  -----------------
13 
14  Project Head: Karl Rupp rupp@iue.tuwien.ac.at
15 
16  (A list of authors and contributors can be found in the PDF manual)
17 
18  License: MIT (X11), see file LICENSE in the base directory
19 ============================================================================= */
20 
25 #include "viennacl/forwards.h"
26 #include "viennacl/scalar.hpp"
27 #include "viennacl/vector.hpp"
28 #include "viennacl/tools/tools.hpp"
30 
31 
32 namespace viennacl
33 {
34  namespace linalg
35  {
36  namespace cuda
37  {
38 
39  namespace detail
40  {
41 
42  template <typename T>
44  const unsigned int * row_index_array,
45  const unsigned int * row_indices,
46  const unsigned int * column_indices,
47  const T * elements,
48  T * vec,
49  unsigned int size)
50  {
51  for (unsigned int row = blockDim.x * blockIdx.x + threadIdx.x;
52  row < size;
53  row += gridDim.x * blockDim.x)
54  {
55  unsigned int eq_row = row_index_array[row];
56  T vec_entry = vec[eq_row];
57  unsigned int row_end = row_indices[row+1];
58 
59  for (unsigned int j = row_indices[row]; j < row_end; ++j)
60  vec_entry -= vec[column_indices[j]] * elements[j];
61 
62  vec[eq_row] = vec_entry;
63  }
64  }
65 
66 
67 
68  template <typename ScalarType>
70  viennacl::backend::mem_handle const & row_index_array,
71  viennacl::backend::mem_handle const & row_buffer,
72  viennacl::backend::mem_handle const & col_buffer,
73  viennacl::backend::mem_handle const & element_buffer,
74  vcl_size_t num_rows
75  )
76  {
77  level_scheduling_substitute_kernel<<<128, 128>>>(detail::cuda_arg<unsigned int>(row_index_array.cuda_handle()),
78  detail::cuda_arg<unsigned int>(row_buffer.cuda_handle()),
79  detail::cuda_arg<unsigned int>(col_buffer.cuda_handle()),
80  detail::cuda_arg<ScalarType>(element_buffer.cuda_handle()),
81  detail::cuda_arg<ScalarType>(vec),
82  static_cast<unsigned int>(num_rows)
83  );
84  }
85 
86  }
87 
88  } // namespace cuda
89  } //namespace linalg
90 } //namespace viennacl
91 
92 
93 #endif
std::size_t vcl_size_t
Definition: forwards.h:58
Common routines for CUDA execution.
Various little tools used here and there in ViennaCL.
__global__ void level_scheduling_substitute_kernel(const unsigned int *row_index_array, const unsigned int *row_indices, const unsigned int *column_indices, const T *elements, T *vec, unsigned int size)
Definition: misc_operations.hpp:43
This file provides the forward declarations for the main types used within ViennaCL.
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:144
void level_scheduling_substitute(vector< ScalarType > &vec, viennacl::backend::mem_handle const &row_index_array, viennacl::backend::mem_handle const &row_buffer, viennacl::backend::mem_handle const &col_buffer, viennacl::backend::mem_handle const &element_buffer, vcl_size_t num_rows)
Definition: misc_operations.hpp:69
vector_expression< const matrix_base< NumericT, F >, const unsigned int, op_row > row(const matrix_base< NumericT, F > &A, unsigned int i)
Definition: matrix.hpp:910
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
Main abstraction class for multiple memory domains. Represents a buffer in either main RAM...
Definition: mem_handle.hpp:62
Implementation of the ViennaCL scalar class.