46 #include <cusp/gallery/poisson.h> 47 #include <cusp/csr_matrix.h> 52 #include "Teuchos_CommandLineProcessor.hpp" 53 #include "Teuchos_TimeMonitor.hpp" 54 #include "Teuchos_StandardCatchMacros.hpp" 56 template <
typename IndexType,
typename ValueType,
typename MemorySpace,
59 IndexType max_its, ValueType tol,
bool verbose) {
61 cusp::csr_matrix<IndexType, ValueType, MemorySpace> A;
64 cusp::gallery::poisson27pt(A, nx, ny, nz);
65 std::cout <<
"N =" << A.num_rows<< std::endl;
66 std::cout <<
"nnz of A = " << A.num_entries << std::endl;
69 std::cout <<
"\nSolving Ax = b with multiple RHS..." << std::endl;
72 cusp::array2d<ValueType, MemorySpace, Orientation>
x(A.num_rows, nrhs, 0);
73 cusp::array2d<ValueType, MemorySpace, Orientation> b(A.num_rows, nrhs, 1);
75 std::cout <<
"numRHS = " << nrhs << std::endl;
85 std::cout <<
"\nPreconditioner statistics" << std::endl;
91 TEUCHOS_FUNC_TIME_MONITOR(
"Total Block-CG Solve Time");
104 typedef int IndexType;
105 typedef double ValueType;
106 typedef cusp::device_memory MemorySpace;
110 bool verbose =
false;
114 Teuchos::CommandLineProcessor
CLP;
116 "This example runs AMG-preconditioned block-CG with CUSP.\n");
119 CLP.setOption(
"nx", &nx,
"Number of mesh points in the x-direction");
121 CLP.setOption(
"ny", &ny,
"Number of mesh points in the y-direction");
123 CLP.setOption(
"nz", &nz,
"Number of mesh points in the z-direction");
125 CLP.setOption(
"nrhs", &nrhs,
"Number of right-hand-sides");
128 "Orientation of block RHS");
129 IndexType max_its = 100;
130 CLP.setOption(
"max_iterations", &max_its,
131 "Maximum number of CG iterations");
133 CLP.setOption(
"tolerance", &tol,
"Convergence tolerance");
135 CLP.setOption(
"device", &device_id,
"CUDA device ID");
136 CLP.setOption(
"verbose",
"quiet", &verbose,
"Verbose output");
140 cudaSetDevice(device_id);
141 cudaDeviceSetSharedMemConfig(cudaSharedMemBankSizeEightByte);
144 cusp_sa_block_cg<IndexType,ValueType,MemorySpace,cusp::row_major>(
145 nx, ny, nz, nrhs, max_its, tol, verbose);
146 else if (orient ==
COL)
147 cusp_sa_block_cg<IndexType,ValueType,MemorySpace,cusp::column_major>(
148 nx, ny, nz, nrhs, max_its, tol, verbose);
150 Teuchos::TimeMonitor::summarize(std::cout);
151 Teuchos::TimeMonitor::zeroOutTimers();
153 TEUCHOS_STANDARD_CATCH_STATEMENTS(verbose, std::cerr, success);
void blockcg(LinearOperator &A, Vector &x, Vector &b)
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
const Orient orient_values[]
int main(int argc, char *argv[])
void cusp_sa_block_cg(IndexType nx, IndexType ny, IndexType nz, IndexType nrhs, IndexType max_its, ValueType tol, bool verbose)
const char * orient_names[]