29 #include <boost/numeric/ublas/matrix_sparse.hpp>
30 #include <boost/numeric/ublas/io.hpp>
31 #include <boost/numeric/ublas/operation_sparse.hpp>
33 #define VIENNACL_WITH_UBLAS 1
53 #ifdef VIENNACL_WITH_OPENCL
65 using namespace boost::numeric;
67 #define BENCHMARK_RUNS 1
70 template<
typename ScalarType>
73 ublas::vector<ScalarType> v2_cpu(v2.
size());
76 for (
unsigned int i=0;i<v1.size(); ++i)
78 if (
std::max( fabs(v2_cpu[i]), fabs(v1[i]) ) > 0 )
79 v2_cpu[i] = fabs(v2_cpu[i] - v1[i]) /
std::max( fabs(v2_cpu[i]), fabs(v1[i]) );
87 template<
typename ScalarType>
90 ublas::vector<ScalarType> v2_cpu(v2.
size());
97 template<
typename MatrixType,
typename VectorType,
typename SolverTag,
typename PrecondTag>
98 void run_solver(MatrixType
const & matrix, VectorType
const & rhs, VectorType
const & ref_result, SolverTag
const & solver, PrecondTag
const & precond,
long ops)
101 VectorType result(rhs);
102 VectorType residual(rhs);
111 double exec_time = timer.
get();
112 std::cout <<
"Exec. time: " << exec_time << std::endl;
113 std::cout <<
"Est. ";
printOps(static_cast<double>(ops), exec_time / BENCHMARK_RUNS);
116 std::cout <<
"Estimated rel. residual: " << solver.error() << std::endl;
117 std::cout <<
"Iterations: " << solver.iters() << std::endl;
118 result -= ref_result;
123 template<
typename ScalarType>
134 ublas::vector<ScalarType> ublas_vec1;
135 ublas::vector<ScalarType> ublas_vec2;
136 ublas::vector<ScalarType> ublas_result;
137 unsigned int solver_iters = 100;
138 unsigned int solver_krylov_dim = 20;
139 double solver_tolerance = 1e-6;
141 ublas::compressed_matrix<ScalarType> ublas_matrix;
144 std::cout <<
"Error reading Matrix file" << std::endl;
148 std::cout <<
"done reading matrix" << std::endl;
150 ublas_result = ublas::scalar_vector<ScalarType>(ublas_matrix.size1(),
ScalarType(1.0));
151 ublas_vec1 =
ublas::prod(ublas_matrix, ublas_result);
152 ublas_vec2 = ublas_vec1;
176 std::cout <<
"------- Jacobi preconditioner ----------" << std::endl;
181 std::cout <<
"------- Row-Scaling preconditioner ----------" << std::endl;
189 std::cout <<
"------- ICHOL0 on CPU (ublas) ----------" << std::endl;
193 exec_time = timer.
get();
194 std::cout <<
"Setup time: " << exec_time << std::endl;
198 ublas_ichol0.
apply(ublas_vec1);
199 exec_time = timer.
get();
200 std::cout <<
"ublas time: " << exec_time << std::endl;
202 std::cout <<
"------- ICHOL0 with ViennaCL ----------" << std::endl;
206 exec_time = timer.
get();
207 std::cout <<
"Setup time: " << exec_time << std::endl;
212 vcl_ichol0.
apply(vcl_vec1);
214 exec_time = timer.
get();
215 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
221 std::cout <<
"------- ILU0 on with ublas ----------" << std::endl;
225 exec_time = timer.
get();
226 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
229 ublas_ilu0.
apply(ublas_vec1);
230 exec_time = timer.
get();
231 std::cout <<
"ublas ILU0 substitution time (no level scheduling): " << exec_time << std::endl;
234 std::cout <<
"------- ILU0 with ViennaCL ----------" << std::endl;
238 exec_time = timer.
get();
239 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
244 vcl_ilu0.
apply(vcl_vec1);
246 exec_time = timer.
get();
247 std::cout <<
"ViennaCL ILU0 substitution time (no level scheduling): " << exec_time << std::endl;
252 exec_time = timer.
get();
253 std::cout <<
"Setup time (with level scheduling): " << exec_time << std::endl;
258 vcl_ilu0_level_scheduling.
apply(vcl_vec1);
260 exec_time = timer.
get();
261 std::cout <<
"ViennaCL ILU0 substitution time (with level scheduling): " << exec_time << std::endl;
267 std::cout <<
"------- Block-ILU0 with ublas ----------" << std::endl;
269 ublas_vec1 = ublas_vec2;
275 exec_time = timer.
get();
276 std::cout <<
"Setup time: " << exec_time << std::endl;
280 ublas_block_ilu0.
apply(ublas_vec1);
281 exec_time = timer.
get();
282 std::cout <<
"ublas time: " << exec_time << std::endl;
284 std::cout <<
"------- Block-ILU0 with ViennaCL ----------" << std::endl;
289 exec_time = timer.
get();
290 std::cout <<
"Setup time: " << exec_time << std::endl;
296 vcl_block_ilu0.
apply(vcl_vec1);
298 exec_time = timer.
get();
299 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
303 std::cout <<
"------- ILUT with ublas ----------" << std::endl;
305 ublas_vec1 = ublas_vec2;
310 exec_time = timer.
get();
311 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
314 ublas_ilut.
apply(ublas_vec1);
315 exec_time = timer.
get();
316 std::cout <<
"ublas ILUT substitution time (no level scheduling): " << exec_time << std::endl;
319 std::cout <<
"------- ILUT with ViennaCL ----------" << std::endl;
323 exec_time = timer.
get();
324 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
329 vcl_ilut.
apply(vcl_vec1);
331 exec_time = timer.
get();
332 std::cout <<
"ViennaCL ILUT substitution time (no level scheduling): " << exec_time << std::endl;
337 exec_time = timer.
get();
338 std::cout <<
"Setup time (with level scheduling): " << exec_time << std::endl;
343 vcl_ilut_level_scheduling.
apply(vcl_vec1);
345 exec_time = timer.
get();
346 std::cout <<
"ViennaCL ILUT substitution time (with level scheduling): " << exec_time << std::endl;
351 std::cout <<
"------- Block-ILUT with ublas ----------" << std::endl;
353 ublas_vec1 = ublas_vec2;
359 exec_time = timer.
get();
360 std::cout <<
"Setup time: " << exec_time << std::endl;
365 ublas_block_ilut.
apply(ublas_vec1);
366 exec_time = timer.
get();
367 std::cout <<
"ublas time: " << exec_time << std::endl;
369 std::cout <<
"------- Block-ILUT with ViennaCL ----------" << std::endl;
374 exec_time = timer.
get();
375 std::cout <<
"Setup time: " << exec_time << std::endl;
381 vcl_block_ilut.
apply(vcl_vec1);
383 exec_time = timer.
get();
384 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
390 long cg_ops =
static_cast<long>(solver_iters * (ublas_matrix.nnz() + 6 * ublas_vec2.size()));
394 std::cout <<
"------- CG solver (no preconditioner) using ublas ----------" << std::endl;
397 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
400 #ifdef VIENNACL_WITH_OPENCL
401 bool is_double = (
sizeof(
ScalarType) ==
sizeof(
double));
404 std::cout <<
"------- CG solver, mixed precision (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
412 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
415 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, ell_matrix ----------" << std::endl;
418 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, sliced_ell_matrix ----------" << std::endl;
421 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, hyb_matrix ----------" << std::endl;
424 std::cout <<
"------- CG solver (ICHOL0 preconditioner) using ublas ----------" << std::endl;
425 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ichol0, cg_ops);
427 std::cout <<
"------- CG solver (ICHOL0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
428 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ichol0, cg_ops);
431 std::cout <<
"------- CG solver (ILU0 preconditioner) using ublas ----------" << std::endl;
432 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ilu0, cg_ops);
434 std::cout <<
"------- CG solver (ILU0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
435 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilu0, cg_ops);
438 std::cout <<
"------- CG solver (Block-ILU0 preconditioner) using ublas ----------" << std::endl;
439 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_block_ilu0, cg_ops);
441 std::cout <<
"------- CG solver (Block-ILU0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
442 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_block_ilu0, cg_ops);
444 std::cout <<
"------- CG solver (ILUT preconditioner) using ublas ----------" << std::endl;
445 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ilut, cg_ops);
447 std::cout <<
"------- CG solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
448 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilut, cg_ops);
450 std::cout <<
"------- CG solver (ILUT preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
451 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilut, cg_ops);
453 std::cout <<
"------- CG solver (Block-ILUT preconditioner) using ublas ----------" << std::endl;
454 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_block_ilut, cg_ops);
456 std::cout <<
"------- CG solver (Block-ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
457 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_block_ilut, cg_ops);
459 std::cout <<
"------- CG solver (Jacobi preconditioner) using ublas ----------" << std::endl;
460 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_jacobi, cg_ops);
462 std::cout <<
"------- CG solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
463 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_jacobi_csr, cg_ops);
465 std::cout <<
"------- CG solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
466 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_jacobi_coo, cg_ops);
469 std::cout <<
"------- CG solver (row scaling preconditioner) using ublas ----------" << std::endl;
470 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_row_scaling, cg_ops);
472 std::cout <<
"------- CG solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
473 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_row_scaling_csr, cg_ops);
475 std::cout <<
"------- CG solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
476 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_row_scaling_coo, cg_ops);
483 long bicgstab_ops =
static_cast<long>(solver_iters * (2 * ublas_matrix.nnz() + 13 * ublas_vec2.size()));
487 std::cout <<
"------- BiCGStab solver (no preconditioner) using ublas ----------" << std::endl;
490 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
493 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
496 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, ell_matrix ----------" << std::endl;
499 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, sliced_ell_matrix ----------" << std::endl;
502 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, hyb_matrix ----------" << std::endl;
505 std::cout <<
"------- BiCGStab solver (ILUT preconditioner) using ublas ----------" << std::endl;
506 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_ilut, bicgstab_ops);
508 std::cout <<
"------- BiCGStab solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
509 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_ilut, bicgstab_ops);
511 std::cout <<
"------- BiCGStab solver (Block-ILUT preconditioner) using ublas ----------" << std::endl;
512 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_block_ilut, bicgstab_ops);
514 #ifdef VIENNACL_WITH_OPENCL
515 std::cout <<
"------- BiCGStab solver (Block-ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
516 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_block_ilut, bicgstab_ops);
522 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) using ublas ----------" << std::endl;
523 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_jacobi, bicgstab_ops);
525 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
526 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_jacobi_csr, bicgstab_ops);
528 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
529 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_jacobi_coo, bicgstab_ops);
532 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) using ublas ----------" << std::endl;
533 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_row_scaling, bicgstab_ops);
535 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
536 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_row_scaling_csr, bicgstab_ops);
538 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
539 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_row_scaling_coo, bicgstab_ops);
546 long gmres_ops =
static_cast<long>(solver_iters * (ublas_matrix.nnz() + (solver_iters * 2 + 7) * ublas_vec2.size()));
550 std::cout <<
"------- GMRES solver (no preconditioner) using ublas ----------" << std::endl;
553 std::cout <<
"------- GMRES solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
556 std::cout <<
"------- GMRES solver (no preconditioner) on GPU, coordinate_matrix ----------" << std::endl;
559 std::cout <<
"------- GMRES solver (no preconditioner) on GPU, ell_matrix ----------" << std::endl;
562 std::cout <<
"------- GMRES solver (no preconditioner) on GPU, sliced_ell_matrix ----------" << std::endl;
565 std::cout <<
"------- GMRES solver (no preconditioner) on GPU, hyb_matrix ----------" << std::endl;
568 std::cout <<
"------- GMRES solver (ILUT preconditioner) using ublas ----------" << std::endl;
569 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_ilut, gmres_ops);
571 std::cout <<
"------- GMRES solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
572 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_ilut, gmres_ops);
574 std::cout <<
"------- GMRES solver (ILUT preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
575 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_ilut, gmres_ops);
578 std::cout <<
"------- GMRES solver (Jacobi preconditioner) using ublas ----------" << std::endl;
579 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_jacobi, gmres_ops);
581 std::cout <<
"------- GMRES solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
582 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_jacobi_csr, gmres_ops);
584 std::cout <<
"------- GMRES solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
585 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_jacobi_coo, gmres_ops);
588 std::cout <<
"------- GMRES solver (row scaling preconditioner) using ublas ----------" << std::endl;
589 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_row_scaling, gmres_ops);
591 std::cout <<
"------- GMRES solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
592 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_row_scaling_csr, gmres_ops);
594 std::cout <<
"------- GMRES solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
595 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_row_scaling_coo, gmres_ops);
602 std::cout << std::endl;
603 std::cout <<
"----------------------------------------------" << std::endl;
604 std::cout <<
" Device Info" << std::endl;
605 std::cout <<
"----------------------------------------------" << std::endl;
607 #ifdef VIENNACL_WITH_OPENCL
609 std::vector<viennacl::ocl::device>
const & devices = pf.
devices();
615 if (devices.size() > 1)
626 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
627 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
628 std::cout <<
" Benchmark for Execution Times of Iterative Solvers provided with ViennaCL " << std::endl;
629 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
630 std::cout <<
" Note that the purpose of this benchmark is not to run solvers until" << std::endl;
631 std::cout <<
" convergence. Instead, only the execution times of a few iterations are" << std::endl;
632 std::cout <<
" recorded. Residual errors are only printed for information." << std::endl << std::endl;
635 std::cout << std::endl;
636 std::cout <<
"----------------------------------------------" << std::endl;
637 std::cout <<
"----------------------------------------------" << std::endl;
638 std::cout <<
"## Benchmark :: Solver" << std::endl;
639 std::cout <<
"----------------------------------------------" << std::endl;
640 std::cout << std::endl;
641 std::cout <<
" -------------------------------" << std::endl;
642 std::cout <<
" # benchmarking single-precision" << std::endl;
643 std::cout <<
" -------------------------------" << std::endl;
644 run_benchmark<float>(ctx);
645 #ifdef VIENNACL_WITH_OPENCL
649 std::cout << std::endl;
650 std::cout <<
" -------------------------------" << std::endl;
651 std::cout <<
" # benchmarking double-precision" << std::endl;
652 std::cout <<
" -------------------------------" << std::endl;
653 run_benchmark<double>(ctx);
Sparse matrix class using a hybrid format composed of the ELL and CSR format for storing the nonzeros...
T norm_2(std::vector< T, A > const &v1)
A reader and writer for the matrix market format is implemented here.
std::vector< platform > get_platforms()
ILU0 preconditioner class, can be supplied to solve()-routines.
This class represents a single scalar value on the GPU and behaves mostly like a built-in scalar type...
Incomplete Cholesky preconditioner class with static pattern (ICHOL0), can be supplied to solve()-rou...
int run_benchmark(viennacl::context ctx)
A tag for incomplete Cholesky factorization with static pattern (ILU0)
Jacobi preconditioner class, can be supplied to solve()-routines. Generic version for non-ViennaCL ma...
A tag for incomplete LU factorization with static pattern (ILU0)
The stabilized bi-conjugate gradient method is implemented here.
ScalarType diff_2(ublas::vector< ScalarType > &v1, viennacl::vector< ScalarType > &v2)
Implementations of incomplete Cholesky factorization preconditioners with static nonzero pattern...
ScalarType diff_inf(ublas::vector< ScalarType > &v1, viennacl::vector< ScalarType > &v2)
void finish()
Synchronizes the execution. finish() will only return after all compute kernels (CUDA, OpenCL) have completed.
bool use_level_scheduling() const
A tag for a jacobi preconditioner.
A block ILU preconditioner class, can be supplied to solve()-routines.
T max(const T &lhs, const T &rhs)
Maximum.
Implementation of a simple Jacobi preconditioner.
Jacobi-type preconditioner class, can be supplied to solve()-routines. This is a diagonal preconditio...
viennacl::ocl::device const & current_device()
Convenience function for returning the active device in the current context.
void printOps(double num_ops, double exec_time)
A tag for the solver GMRES. Used for supplying solver parameters and for dispatching the solve() func...
Implementation of the coordinate_matrix class.
Implementation of a OpenCL-like context, which serves as a unification of {OpenMP, CUDA, OpenCL} at the user API.
std::string info(vcl_size_t indent=0, char indent_char= ' ') const
Returns an info string with a few properties of the device. Use full_info() to get all details...
Represents a generic 'context' similar to an OpenCL context, but is backend-agnostic and thus also su...
viennacl::vector< NumericT > solve(viennacl::compressed_matrix< NumericT > const &A, viennacl::vector_base< NumericT > const &rhs, bicgstab_tag const &tag, viennacl::linalg::no_precond)
Overload for the pipelined BiCGStab implementation for the ViennaCL sparse matrix types...
void apply(VectorT &vec) const
viennacl::vector< float > v1
Implementation of the hyb_matrix class.
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
Implementations of the generalized minimum residual method are in this file.
Sparse matrix class using the ELLPACK format for storing the nonzeros.
iterator begin()
Returns an iterator pointing to the beginning of the vector (STL like)
A tag class representing the use of no preconditioner.
Implementations of incomplete factorization preconditioners. Convenience header file.
A tag for incomplete LU factorization with threshold (ILUT)
Sparse matrix class using the sliced ELLPACK with parameters C, .
Implementation of the compressed_matrix class.
Implementation of the sliced_ell_matrix class.
A tag for a row scaling preconditioner which merely normalizes the equation system such that each row...
void run_solver(MatrixType const &matrix, VectorType const &rhs, VectorType const &ref_result, SolverTag const &solver, PrecondTag const &precond, long ops)
bool double_support() const
ViennaCL convenience function: Returns true if the device supports double precision.
ILUT preconditioner class, can be supplied to solve()-routines.
The conjugate gradient method is implemented here.
void apply(VectorT &vec) const
Implementation of the ell_matrix class.
void prod(const MatrixT1 &A, bool transposed_A, const MatrixT2 &B, bool transposed_B, MatrixT3 &C, ScalarT alpha, ScalarT beta)
void apply(VectorT &vec) const
A row normalization preconditioner is implemented here.
void apply(VectorT &vec) const
viennacl::vector< int > v2
bool use_level_scheduling() const
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
T norm_inf(std::vector< T, A > const &v1)
void copy(std::vector< NumericT > &cpu_vec, circulant_matrix< NumericT, AlignmentV > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
size_type size() const
Returns the length of the vector (cf. std::vector)
A tag for the conjugate gradient Used for supplying solver parameters and for dispatching the solve()...
A sparse square matrix in compressed sparse rows format.
A tag for the stabilized Bi-conjugate gradient solver. Used for supplying solver parameters and for d...
long read_matrix_market_file(MatrixT &mat, const char *file, long index_base=1)
Reads a sparse matrix from a file (MatrixMarket format)
iterator end()
Returns an iterator pointing to the end of the vector (STL like)
viennacl::ocl::context & get_context(long i)
Convenience function for returning the current context.
A tag for the conjugate gradient Used for supplying solver parameters and for dispatching the solve()...
Implementation of the ViennaCL scalar class.
void setup_context(long i, std::vector< cl_device_id > const &devices)
Convenience function for setting devices for a context.
The conjugate gradient method using mixed precision is implemented here. Experimental.
A sparse square matrix, where entries are stored as triplets (i,j, val), where i and j are the row an...