ViennaCL - The Vienna Computing Library  1.6.2
Free open-source GPU-accelerated linear algebra and solver library.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
sliced_ell_matrix.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_OPENCL_KERNELS_SLICED_ELL_MATRIX_HPP
2 #define VIENNACL_LINALG_OPENCL_KERNELS_SLICED_ELL_MATRIX_HPP
3 
7 #include "viennacl/ocl/utils.hpp"
8 
10 
13 namespace viennacl
14 {
15 namespace linalg
16 {
17 namespace opencl
18 {
19 namespace kernels
20 {
21 
23 
24 template<typename StringT>
25 void generate_sliced_ell_vec_mul(StringT & source, std::string const & numeric_string)
26 {
27  source.append("__kernel void vec_mul( \n");
28  source.append(" __global const unsigned int * columns_per_block, \n");
29  source.append(" __global const unsigned int * column_indices, \n");
30  source.append(" __global const unsigned int * block_start, \n");
31  source.append(" __global const "); source.append(numeric_string); source.append(" * elements, \n");
32  source.append(" __global const "); source.append(numeric_string); source.append(" * x, \n");
33  source.append(" uint4 layout_x, \n");
34  source.append(" __global "); source.append(numeric_string); source.append(" * result, \n");
35  source.append(" uint4 layout_result) \n");
36  source.append("{ \n");
37  source.append(" uint local_id = get_local_id(0); \n");
38  source.append(" uint local_size = get_local_size(0); \n");
39  source.append(" uint num_rows = layout_result.z; \n");
40 
41  source.append(" for (uint block_idx = get_group_id(0); block_idx <= num_rows / local_size; block_idx += get_num_groups(0)) { \n");
42  source.append(" "); source.append(numeric_string); source.append(" sum = 0; \n");
43 
44  source.append(" uint row = block_idx * local_size + local_id; \n");
45  source.append(" uint offset = block_start[block_idx]; \n");
46  source.append(" uint num_columns = columns_per_block[block_idx]; \n");
47  source.append(" for (uint item_id = 0; item_id < num_columns; item_id++) { \n");
48  source.append(" uint index = offset + item_id * local_size + local_id; \n");
49  source.append(" "); source.append(numeric_string); source.append(" val = elements[index]; \n");
50  source.append(" sum += val ? (x[column_indices[index] * layout_x.y + layout_x.x] * val) : 0; \n");
51  source.append(" } \n");
52 
53  source.append(" if (row < num_rows) \n");
54  source.append(" result[row * layout_result.y + layout_result.x] = sum; \n");
55  source.append(" } \n");
56  source.append("} \n");
57 }
58 
59 
61 
62 // main kernel class
64 template<typename NumericT, typename IndexT>
66 
67 template<typename NumericT>
68 struct sliced_ell_matrix<NumericT, unsigned int>
69 {
70  static std::string program_name()
71  {
73  }
74 
75  static void init(viennacl::ocl::context & ctx)
76  {
77  static std::map<cl_context, bool> init_done;
78  if (!init_done[ctx.handle().get()])
79  {
81  std::string numeric_string = viennacl::ocl::type_to_string<NumericT>::apply();
82 
83  std::string source;
84  source.reserve(1024);
85 
86  viennacl::ocl::append_double_precision_pragma<NumericT>(ctx, source);
87 
88  // fully parametrized kernels:
89  generate_sliced_ell_vec_mul(source, numeric_string);
90 
91  std::string prog_name = program_name();
92  #ifdef VIENNACL_BUILD_INFO
93  std::cout << "Creating program " << prog_name << std::endl;
94  #endif
95  ctx.add_program(source, prog_name);
96  init_done[ctx.handle().get()] = true;
97  } //if
98  } //init
99 };
100 
101 } // namespace kernels
102 } // namespace opencl
103 } // namespace linalg
104 } // namespace viennacl
105 #endif
106 
Implements a OpenCL platform within ViennaCL.
void generate_sliced_ell_vec_mul(StringT &source, std::string const &numeric_string)
Various little tools used here and there in ViennaCL.
Manages an OpenCL context and provides the respective convenience functions for creating buffers...
Definition: context.hpp:54
Provides OpenCL-related utilities.
const viennacl::ocl::handle< cl_context > & handle() const
Returns the context handle.
Definition: context.hpp:613
Common implementations shared by OpenCL-based operations.
static void apply(viennacl::ocl::context const &)
Definition: utils.hpp:40
const OCL_TYPE & get() const
Definition: handle.hpp:189
Main kernel class for generating OpenCL kernels for ell_matrix.
Representation of an OpenCL kernel in ViennaCL.
Helper class for converting a type to its string representation.
Definition: utils.hpp:57