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
52 #ifdef VIENNACL_WITH_OPENCL
66 #define BENCHMARK_RUNS 1
69 template<
typename ScalarType>
72 ublas::vector<ScalarType> v2_cpu(v2.
size());
75 for (
unsigned int i=0;i<v1.size(); ++i)
77 if (
std::max( fabs(v2_cpu[i]), fabs(v1[i]) ) > 0 )
78 v2_cpu[i] = fabs(v2_cpu[i] - v1[i]) /
std::max( fabs(v2_cpu[i]), fabs(v1[i]) );
86 template<
typename ScalarType>
89 ublas::vector<ScalarType> v2_cpu(v2.
size());
96 template<
typename MatrixType,
typename VectorType,
typename SolverTag,
typename PrecondTag>
97 void run_solver(MatrixType
const & matrix, VectorType
const & rhs, VectorType
const & ref_result, SolverTag
const & solver, PrecondTag
const & precond,
long ops)
100 VectorType result(rhs);
101 VectorType residual(rhs);
110 double exec_time = timer.
get();
111 std::cout <<
"Exec. time: " << exec_time << std::endl;
112 std::cout <<
"Est. ";
printOps(static_cast<double>(ops), exec_time / BENCHMARK_RUNS);
115 std::cout <<
"Estimated rel. residual: " << solver.error() << std::endl;
116 std::cout <<
"Iterations: " << solver.iters() << std::endl;
117 result -= ref_result;
122 template<
typename ScalarType>
133 ublas::vector<ScalarType> ublas_vec1;
134 ublas::vector<ScalarType> ublas_vec2;
135 ublas::vector<ScalarType> ublas_result;
136 unsigned int solver_iters = 100;
137 unsigned int solver_krylov_dim = 20;
138 double solver_tolerance = 1e-6;
140 ublas::compressed_matrix<ScalarType> ublas_matrix;
143 std::cout <<
"Error reading Matrix file" << std::endl;
147 std::cout <<
"done reading matrix" << std::endl;
149 ublas_result = ublas::scalar_vector<ScalarType>(ublas_matrix.size1(),
ScalarType(1.0));
150 ublas_vec1 =
ublas::prod(ublas_matrix, ublas_result);
151 ublas_vec2 = ublas_vec1;
173 std::cout <<
"------- Jacobi preconditioner ----------" << std::endl;
178 std::cout <<
"------- Row-Scaling preconditioner ----------" << std::endl;
186 std::cout <<
"------- ICHOL0 on CPU (ublas) ----------" << std::endl;
190 exec_time = timer.
get();
191 std::cout <<
"Setup time: " << exec_time << std::endl;
195 ublas_ichol0.
apply(ublas_vec1);
196 exec_time = timer.
get();
197 std::cout <<
"ublas time: " << exec_time << std::endl;
199 std::cout <<
"------- ICHOL0 with ViennaCL ----------" << std::endl;
203 exec_time = timer.
get();
204 std::cout <<
"Setup time: " << exec_time << std::endl;
209 vcl_ichol0.
apply(vcl_vec1);
211 exec_time = timer.
get();
212 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
218 std::cout <<
"------- ILU0 on with ublas ----------" << std::endl;
222 exec_time = timer.
get();
223 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
226 ublas_ilu0.
apply(ublas_vec1);
227 exec_time = timer.
get();
228 std::cout <<
"ublas ILU0 substitution time (no level scheduling): " << exec_time << std::endl;
231 std::cout <<
"------- ILU0 with ViennaCL ----------" << std::endl;
235 exec_time = timer.
get();
236 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
241 vcl_ilu0.
apply(vcl_vec1);
243 exec_time = timer.
get();
244 std::cout <<
"ViennaCL ILU0 substitution time (no level scheduling): " << exec_time << std::endl;
249 exec_time = timer.
get();
250 std::cout <<
"Setup time (with level scheduling): " << exec_time << std::endl;
255 vcl_ilu0_level_scheduling.
apply(vcl_vec1);
257 exec_time = timer.
get();
258 std::cout <<
"ViennaCL ILU0 substitution time (with level scheduling): " << exec_time << std::endl;
264 std::cout <<
"------- Block-ILU0 with ublas ----------" << std::endl;
266 ublas_vec1 = ublas_vec2;
272 exec_time = timer.
get();
273 std::cout <<
"Setup time: " << exec_time << std::endl;
277 ublas_block_ilu0.
apply(ublas_vec1);
278 exec_time = timer.
get();
279 std::cout <<
"ublas time: " << exec_time << std::endl;
281 std::cout <<
"------- Block-ILU0 with ViennaCL ----------" << std::endl;
286 exec_time = timer.
get();
287 std::cout <<
"Setup time: " << exec_time << std::endl;
293 vcl_block_ilu0.
apply(vcl_vec1);
295 exec_time = timer.
get();
296 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
300 std::cout <<
"------- ILUT with ublas ----------" << std::endl;
302 ublas_vec1 = ublas_vec2;
307 exec_time = timer.
get();
308 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
311 ublas_ilut.
apply(ublas_vec1);
312 exec_time = timer.
get();
313 std::cout <<
"ublas ILUT substitution time (no level scheduling): " << exec_time << std::endl;
316 std::cout <<
"------- ILUT with ViennaCL ----------" << std::endl;
320 exec_time = timer.
get();
321 std::cout <<
"Setup time (no level scheduling): " << exec_time << std::endl;
326 vcl_ilut.
apply(vcl_vec1);
328 exec_time = timer.
get();
329 std::cout <<
"ViennaCL ILUT substitution time (no level scheduling): " << exec_time << std::endl;
334 exec_time = timer.
get();
335 std::cout <<
"Setup time (with level scheduling): " << exec_time << std::endl;
340 vcl_ilut_level_scheduling.
apply(vcl_vec1);
342 exec_time = timer.
get();
343 std::cout <<
"ViennaCL ILUT substitution time (with level scheduling): " << exec_time << std::endl;
348 std::cout <<
"------- Block-ILUT with ublas ----------" << std::endl;
350 ublas_vec1 = ublas_vec2;
356 exec_time = timer.
get();
357 std::cout <<
"Setup time: " << exec_time << std::endl;
362 ublas_block_ilut.
apply(ublas_vec1);
363 exec_time = timer.
get();
364 std::cout <<
"ublas time: " << exec_time << std::endl;
366 std::cout <<
"------- Block-ILUT with ViennaCL ----------" << std::endl;
371 exec_time = timer.
get();
372 std::cout <<
"Setup time: " << exec_time << std::endl;
378 vcl_block_ilut.
apply(vcl_vec1);
380 exec_time = timer.
get();
381 std::cout <<
"ViennaCL time: " << exec_time << std::endl;
387 long cg_ops =
static_cast<long>(solver_iters * (ublas_matrix.nnz() + 6 * ublas_vec2.size()));
391 std::cout <<
"------- CG solver (no preconditioner) using ublas ----------" << std::endl;
394 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
397 #ifdef VIENNACL_WITH_OPENCL
398 bool is_double = (
sizeof(
ScalarType) ==
sizeof(
double));
401 std::cout <<
"------- CG solver, mixed precision (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
409 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
412 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, ell_matrix ----------" << std::endl;
415 std::cout <<
"------- CG solver (no preconditioner) via ViennaCL, hyb_matrix ----------" << std::endl;
418 std::cout <<
"------- CG solver (ICHOL0 preconditioner) using ublas ----------" << std::endl;
419 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ichol0, cg_ops);
421 std::cout <<
"------- CG solver (ICHOL0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
422 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ichol0, cg_ops);
425 std::cout <<
"------- CG solver (ILU0 preconditioner) using ublas ----------" << std::endl;
426 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ilu0, cg_ops);
428 std::cout <<
"------- CG solver (ILU0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
429 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilu0, cg_ops);
432 std::cout <<
"------- CG solver (Block-ILU0 preconditioner) using ublas ----------" << std::endl;
433 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_block_ilu0, cg_ops);
435 std::cout <<
"------- CG solver (Block-ILU0 preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
436 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_block_ilu0, cg_ops);
438 std::cout <<
"------- CG solver (ILUT preconditioner) using ublas ----------" << std::endl;
439 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_ilut, cg_ops);
441 std::cout <<
"------- CG solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
442 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilut, cg_ops);
444 std::cout <<
"------- CG solver (ILUT preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
445 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_ilut, cg_ops);
447 std::cout <<
"------- CG solver (Block-ILUT preconditioner) using ublas ----------" << std::endl;
448 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_block_ilut, cg_ops);
450 std::cout <<
"------- CG solver (Block-ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
451 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_block_ilut, cg_ops);
453 std::cout <<
"------- CG solver (Jacobi preconditioner) using ublas ----------" << std::endl;
454 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_jacobi, cg_ops);
456 std::cout <<
"------- CG solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
457 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_jacobi_csr, cg_ops);
459 std::cout <<
"------- CG solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
460 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_jacobi_coo, cg_ops);
463 std::cout <<
"------- CG solver (row scaling preconditioner) using ublas ----------" << std::endl;
464 run_solver(ublas_matrix, ublas_vec2, ublas_result, cg_solver, ublas_row_scaling, cg_ops);
466 std::cout <<
"------- CG solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
467 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, cg_solver, vcl_row_scaling_csr, cg_ops);
469 std::cout <<
"------- CG solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
470 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, cg_solver, vcl_row_scaling_coo, cg_ops);
477 long bicgstab_ops =
static_cast<long>(solver_iters * (2 * ublas_matrix.nnz() + 13 * ublas_vec2.size()));
481 std::cout <<
"------- BiCGStab solver (no preconditioner) using ublas ----------" << std::endl;
484 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
487 std::cout <<
"------- BiCGStab solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
491 std::cout <<
"------- BiCGStab solver (ILUT preconditioner) using ublas ----------" << std::endl;
492 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_ilut, bicgstab_ops);
494 std::cout <<
"------- BiCGStab solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
495 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_ilut, bicgstab_ops);
497 std::cout <<
"------- BiCGStab solver (Block-ILUT preconditioner) using ublas ----------" << std::endl;
498 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_block_ilut, bicgstab_ops);
500 #ifdef VIENNACL_WITH_OPENCL
501 std::cout <<
"------- BiCGStab solver (Block-ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
502 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_block_ilut, bicgstab_ops);
508 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) using ublas ----------" << std::endl;
509 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_jacobi, bicgstab_ops);
511 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
512 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_jacobi_csr, bicgstab_ops);
514 std::cout <<
"------- BiCGStab solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
515 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_jacobi_coo, bicgstab_ops);
518 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) using ublas ----------" << std::endl;
519 run_solver(ublas_matrix, ublas_vec2, ublas_result, bicgstab_solver, ublas_row_scaling, bicgstab_ops);
521 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
522 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_row_scaling_csr, bicgstab_ops);
524 std::cout <<
"------- BiCGStab solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
525 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, bicgstab_solver, vcl_row_scaling_coo, bicgstab_ops);
532 long gmres_ops =
static_cast<long>(solver_iters * (ublas_matrix.nnz() + (solver_iters * 2 + 7) * ublas_vec2.size()));
536 std::cout <<
"------- GMRES solver (no preconditioner) using ublas ----------" << std::endl;
539 std::cout <<
"------- GMRES solver (no preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
542 std::cout <<
"------- GMRES solver (no preconditioner) on GPU, coordinate_matrix ----------" << std::endl;
546 std::cout <<
"------- GMRES solver (ILUT preconditioner) using ublas ----------" << std::endl;
547 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_ilut, gmres_ops);
549 std::cout <<
"------- GMRES solver (ILUT preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
550 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_ilut, gmres_ops);
552 std::cout <<
"------- GMRES solver (ILUT preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
553 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_ilut, gmres_ops);
556 std::cout <<
"------- GMRES solver (Jacobi preconditioner) using ublas ----------" << std::endl;
557 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_jacobi, gmres_ops);
559 std::cout <<
"------- GMRES solver (Jacobi preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
560 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_jacobi_csr, gmres_ops);
562 std::cout <<
"------- GMRES solver (Jacobi preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
563 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_jacobi_coo, gmres_ops);
566 std::cout <<
"------- GMRES solver (row scaling preconditioner) using ublas ----------" << std::endl;
567 run_solver(ublas_matrix, ublas_vec2, ublas_result, gmres_solver, ublas_row_scaling, gmres_ops);
569 std::cout <<
"------- GMRES solver (row scaling preconditioner) via ViennaCL, compressed_matrix ----------" << std::endl;
570 run_solver(vcl_compressed_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_row_scaling_csr, gmres_ops);
572 std::cout <<
"------- GMRES solver (row scaling preconditioner) via ViennaCL, coordinate_matrix ----------" << std::endl;
573 run_solver(vcl_coordinate_matrix, vcl_vec2, vcl_result, gmres_solver, vcl_row_scaling_coo, gmres_ops);
580 std::cout << std::endl;
581 std::cout <<
"----------------------------------------------" << std::endl;
582 std::cout <<
" Device Info" << std::endl;
583 std::cout <<
"----------------------------------------------" << std::endl;
585 #ifdef VIENNACL_WITH_OPENCL
587 std::vector<viennacl::ocl::device>
const & devices = pf.
devices();
593 if (devices.size() > 1)
604 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
605 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
606 std::cout <<
" Benchmark for Execution Times of Iterative Solvers provided with ViennaCL " << std::endl;
607 std::cout <<
"---------------------------------------------------------------------------" << std::endl;
608 std::cout <<
" Note that the purpose of this benchmark is not to run solvers until" << std::endl;
609 std::cout <<
" convergence. Instead, only the execution times of a few iterations are" << std::endl;
610 std::cout <<
" recorded. Residual errors are only printed for information." << std::endl << std::endl;
613 std::cout << std::endl;
614 std::cout <<
"----------------------------------------------" << std::endl;
615 std::cout <<
"----------------------------------------------" << std::endl;
616 std::cout <<
"## Benchmark :: Solver" << std::endl;
617 std::cout <<
"----------------------------------------------" << std::endl;
618 std::cout << std::endl;
619 std::cout <<
" -------------------------------" << std::endl;
620 std::cout <<
" # benchmarking single-precision" << std::endl;
621 std::cout <<
" -------------------------------" << std::endl;
622 run_benchmark<float>(ctx);
623 #ifdef VIENNACL_WITH_OPENCL
627 std::cout << std::endl;
628 std::cout <<
" -------------------------------" << std::endl;
629 std::cout <<
" # benchmarking double-precision" << std::endl;
630 std::cout <<
" -------------------------------" << std::endl;
631 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...
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)
Implementation of the compressed_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...
viennacl::vector< NumericT > solve(MatrixT const &A, viennacl::vector_base< NumericT > const &rhs, bicgstab_tag const &tag, viennacl::linalg::no_precond)
Implementation of a pipelined stabilized Bi-conjugate gradient solver.