1 #ifndef VIENNACL_LINALG_CUDA_BISECT_KERNEL_CALLS_HPP_
2 #define VIENNACL_LINALG_CUDA_BISECT_KERNEL_CALLS_HPP_
42 template<
typename NumericT>
44 const unsigned int mat_size,
45 const NumericT lg,
const NumericT ug,
46 const NumericT precision)
53 bisectKernelSmall<<< blocks, threads >>>(
54 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_a),
55 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_b) + 1,
57 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
vcl_g_left),
58 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
vcl_g_right),
59 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
vcl_g_left_count),
68 template<
typename NumericT>
70 const unsigned int mat_size,
71 const NumericT lg,
const NumericT ug,
72 const NumericT precision)
77 bisectKernelLarge<<< blocks, threads >>>
78 (viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_a),
79 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_b) + 1,
81 lg, ug, 0, mat_size, precision,
82 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_num_one),
84 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_left_one),
85 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_right_one),
86 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_pos_one),
87 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_left_mult),
88 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_right_mult),
91 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_blocks_mult),
100 template<
typename NumericT>
102 const unsigned int mat_size,
103 const NumericT precision)
106 unsigned int num_one_intervals = result.
g_num_one;
109 grid_onei.x = num_blocks;
110 grid_onei.y = 1, grid_onei.z = 1;
114 bisectKernelLarge_OneIntervals<<< grid_onei , threads_onei >>>
115 (viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_a),
116 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_b) + 1,
117 mat_size, num_one_intervals,
118 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_left_one),
119 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_right_one),
120 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_pos_one),
129 template<
typename NumericT>
131 const unsigned int mat_size,
132 const NumericT precision)
140 dim3 grid_mult(num_blocks_mult, 1, 1);
143 bisectKernelLarge_MultIntervals<<< grid_mult, threads_mult >>>
144 (viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_a),
145 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(input.
g_b) + 1,
147 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_blocks_mult),
148 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_blocks_mult_sum),
149 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_left_mult),
150 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_right_mult),
151 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_left_count_mult),
153 viennacl::linalg::cuda::detail::cuda_arg<NumericT>(result.
g_lambda_mult),
154 viennacl::linalg::cuda::detail::cuda_arg<unsigned int>(result.
g_pos_mult),
void bisectLarge(const viennacl::linalg::detail::InputData< NumericT > &input, viennacl::linalg::detail::ResultDataLarge< NumericT > &result, const unsigned int mat_size, const NumericT lg, const NumericT ug, const NumericT precision)
viennacl::vector< NumericT > g_left_mult
#define VIENNACL_BISECT_MAX_THREADS_BLOCK
Determine eigenvalues for small symmetric, tridiagonal matrix.
viennacl::vector< NumericT > g_left_one
viennacl::scalar< unsigned int > g_num_blocks_mult
viennacl::vector< unsigned int > g_blocks_mult
viennacl::vector< unsigned int > g_left_count_mult
Main namespace in ViennaCL. Holds all the basic types such as vector, matrix, etc. and defines operations upon them.
viennacl::vector< unsigned int > g_right_count_mult
unsigned int getNumBlocksLinear(const unsigned int num_threads, const unsigned int num_threads_block)
viennacl::vector< NumericT > vcl_g_right
right interval limits at the end of the computation
viennacl::vector< NumericT > g_lambda_mult
viennacl::vector< NumericT > vcl_g_left
left interval limits at the end of the computation
viennacl::vector< unsigned int > vcl_g_left_count
number of eigenvalues smaller than the left interval limit
viennacl::vector< unsigned int > g_pos_one
#define VIENNACL_BISECT_MAX_THREADS_BLOCK_SMALL_MATRIX
First step of the bisection algorithm for the computation of eigenvalues.
Determine eigenvalues for large matrices for intervals that contained after the first step one eigenv...
viennacl::vector< NumericT > g_right_one
viennacl::vector< unsigned int > g_pos_mult
void bisectLarge_OneIntervals(const viennacl::linalg::detail::InputData< NumericT > &input, viennacl::linalg::detail::ResultDataLarge< NumericT > &result, const unsigned int mat_size, const NumericT precision)
viennacl::vector< unsigned int > vcl_g_right_count
number of eigenvalues bigger than the right interval limit
In this class the data of the result for large matrices is stored.
void bisectSmall(const viennacl::linalg::detail::InputData< NumericT > &input, viennacl::linalg::detail::ResultDataSmall< NumericT > &result, const unsigned int mat_size, const NumericT lg, const NumericT ug, const NumericT precision)
viennacl::vector< unsigned int > g_blocks_mult_sum
viennacl::scalar< unsigned int > g_num_one
number of intervals containing one eigenvalue after the first step
#define VIENNACL_CUDA_LAST_ERROR_CHECK(message)
void bisectLarge_MultIntervals(const viennacl::linalg::detail::InputData< NumericT > &input, viennacl::linalg::detail::ResultDataLarge< NumericT > &result, const unsigned int mat_size, const NumericT precision)
In this class the data of the result for small matrices is stored.
Second step of the bisection algorithm for the computation of eigenvalues for large matrices...
viennacl::vector< NumericT > g_right_mult