1 #ifndef VIENNACL_LINALG_CUDA_MISC_OPERATIONS_HPP_
2 #define VIENNACL_LINALG_CUDA_MISC_OPERATIONS_HPP_
44 const unsigned int * row_index_array,
45 const unsigned int * row_indices,
46 const unsigned int * column_indices,
51 for (
unsigned int row = blockDim.x * blockIdx.x + threadIdx.x;
53 row += gridDim.x * blockDim.x)
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];
59 for (
unsigned int j = row_indices[
row]; j < row_end; ++j)
60 vec_entry -= vec[column_indices[j]] * elements[j];
62 vec[eq_row] = vec_entry;
68 template <
typename ScalarType>
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)
std::size_t vcl_size_t
Definition: forwards.h:58
Common routines for CUDA execution.
__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.