57 template<
typename MatrixType,
typename VectorType,
typename SolverTag,
typename PrecondTag>
58 void run_solver(MatrixType
const & matrix, VectorType
const & rhs, VectorType
const & ref_result, SolverTag
const & solver, PrecondTag
const & precond)
60 VectorType result(rhs);
61 VectorType residual(rhs);
67 std::cout <<
" > Solver time: " << timer.
get() << std::endl;
70 std::cout <<
" > Iterations: " << solver.iters() << std::endl;
80 template<
typename ScalarType>
88 std::cout <<
"-- CG with AMG preconditioner, " << info <<
" --" << std::endl;
91 std::cout <<
" * Setup phase (ViennaCL types)..." << std::endl;
96 std::cout <<
" > Setup time: " << timer.
get() << std::endl;
98 std::cout <<
" * CG solver (ViennaCL types)..." << std::endl;
99 run_solver(vcl_compressed_matrix, vcl_vec, vcl_result, cg_solver, vcl_amg);
107 int main(
int argc,
char **argv)
109 std::string filename(
"../examples/testdata/mat65k.mtx");
116 std::cout << std::endl;
117 std::cout <<
"----------------------------------------------" << std::endl;
118 std::cout <<
" Device Info" << std::endl;
119 std::cout <<
"----------------------------------------------" << std::endl;
121 #ifdef VIENNACL_WITH_OPENCL
124 std::vector<viennacl::ocl::device>
const & devices = pf.
devices();
130 if (devices.size() > 1)
151 std::cout <<
"Reading matrix..." << std::endl;
152 std::vector< std::map<unsigned int, ScalarType> > read_in_matrix;
155 std::cout <<
"Error reading Matrix file" << std::endl;
159 std::cout <<
"Reading matrix completed." << std::endl;
164 std::vector<ScalarType> std_vec, std_result;
170 for (std::size_t i=0; i<std_result.
size(); ++i)
192 std::cout <<
"-- CG solver (no preconditioner, warmup) --" << std::endl;
207 run_amg(cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"ONEPASS COARSENING, DIRECT INTERPOLATION", amg_tag_direct);
215 run_amg(cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"AG COARSENING (PMIS), AG INTERPOLATION", amg_tag_agg_pmis);
223 run_amg (cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"AG COARSENING (PMIS), SA INTERPOLATION", amg_tag_sa_pmis);
225 std::cout << std::endl;
226 std::cout <<
" -------------- Benchmark runs -------------- " << std::endl;
227 std::cout << std::endl;
229 std::cout <<
"-- CG solver (no preconditioner) --" << std::endl;
231 run_amg(cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"ONEPASS COARSENING, DIRECT INTERPOLATION", amg_tag_direct);
232 run_amg(cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"AG COARSENING (PMIS), AG INTERPOLATION", amg_tag_agg_pmis);
233 run_amg (cg_solver, vcl_vec, vcl_result, vcl_compressed_matrix,
"AG COARSENING (PMIS), SA INTERPOLATION", amg_tag_sa_pmis);
238 std::cout <<
"!!!! TUTORIAL COMPLETED SUCCESSFULLY !!!!" << std::endl;
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()
Generic interface for the l^2-norm. See viennacl/linalg/vector_operations.hpp for implementations...
void set_setup_context(viennacl::context ctx)
Sets the ViennaCL context for the setup stage. Set this to a host context if you want to run the setu...
void set_strong_connection_threshold(double threshold)
Sets the strong connection threshold. Customizations by the user essential for best results! ...
void set_jacobi_weight(double w)
Sets the weight (damping) for the Jacobi smoother.
const vcl_size_t & size1() const
Returns the number of rows.
Helper routines for generating sparse matrices.
void set_target_context(viennacl::context ctx)
Sets the ViennaCL context for the solver cycle stage (i.e. preconditioner applications).
The stabilized bi-conjugate gradient method is implemented here.
void finish()
Synchronizes the execution. finish() will only return after all compute kernels (CUDA, OpenCL) have completed.
void set_postsmooth_steps(vcl_size_t steps)
Sets the number of smoother applications on the coarse level before interpolation to the finer level...
viennacl::ocl::device const & current_device()
Convenience function for returning the active device in the current context.
Implementation of the coordinate_matrix class.
VectorT solve(MatrixT const &matrix, VectorT const &rhs, bicgstab_tag const &tag, PreconditionerT const &precond)
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...
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
void set_interpolation_method(amg_interpolation_method interpol)
Sets the interpolation method to the provided method.
A tag class representing the use of no preconditioner.
Implementations of incomplete factorization preconditioners. Convenience header file.
void set_presmooth_steps(vcl_size_t steps)
Sets the number of smoother applications on the fine level before restriction to the coarser level...
Implementation of the compressed_matrix class.
void run_solver(MatrixType const &matrix, VectorType const &rhs, VectorType const &ref_result, SolverTag const &solver, PrecondTag const &precond, long ops)
Main include file for algebraic multigrid (AMG) preconditioners. Experimental.
AMG preconditioner class, can be supplied to solve()-routines.
The conjugate gradient method is implemented here.
void resize(size_type new_size, bool preserve=true)
Resizes the allocated memory for the vector. Pads the memory to be a multiple of 'AlignmentV'.
void set_coarsening_method(amg_coarsening_method s)
Sets the strategy used for constructing coarse grids.
A simple, yet (mostly) sufficiently accurate timer for benchmarking and profiling.
A tag for algebraic multigrid (AMG). Used to transport information from the user to the implementatio...
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object.
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
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.
long read_matrix_market_file(MatrixT &mat, const char *file, long index_base=1)
Reads a sparse matrix from a file (MatrixMarket format)
viennacl::ocl::context & get_context(long i)
Convenience function for returning the current context.
void setup_context(long i, std::vector< cl_device_id > const &devices)
Convenience function for setting devices for a context.