62 #ifdef HAVE_BELOS_TRIUTILS
63 #include "Trilinos_Util_iohb.h"
70 using namespace Teuchos;
72 int main(
int argc,
char *argv[]) {
75 typedef std::complex<double> ST;
77 typedef std::complex<double> ST;
79 std::cout <<
"Not compiled with std::complex support." << std::endl;
80 std::cout <<
"End Result: TEST FAILED" << std::endl;
85 typedef SCT::magnitudeType MT;
91 ST zero = SCT::zero();
100 bool verbose =
false;
103 bool norm_failure =
false;
104 bool proc_verbose =
false;
105 bool userandomrhs =
true;
111 int maxsubspace = 50;
112 int maxrestarts = 15;
113 std::string outersolver(
"Block Gmres");
114 std::string
filename(
"mhd1280b.cua");
119 cmdp.
setOption(
"verbose",
"quiet",&verbose,
"Print messages and results.");
120 cmdp.
setOption(
"use-random-rhs",
"use-rhs",&userandomrhs,
"Use linear system RHS or random RHS to generate polynomial.");
121 cmdp.
setOption(
"frequency",&frequency,
"Solvers frequency for printing residuals (#iters).");
122 cmdp.
setOption(
"filename",&filename,
"Filename for test matrix. Acceptable file extensions: *.hb,*.mtx,*.triU,*.triS");
123 cmdp.
setOption(
"outersolver",&outersolver,
"Name of outer solver to be used with GMRES poly");
124 cmdp.
setOption(
"tol",&tol,
"Relative residual tolerance used by GMRES solver.");
125 cmdp.
setOption(
"poly-tol",&polytol,
"Relative residual tolerance used to construct the GMRES polynomial.");
126 cmdp.
setOption(
"num-rhs",&numrhs,
"Number of right-hand sides to be solved for.");
127 cmdp.
setOption(
"block-size",&blocksize,
"Block size used by GMRES.");
128 cmdp.
setOption(
"max-iters",&maxiters,
"Maximum number of iterations per linear system (-1 = adapted to problem/block size).");
129 cmdp.
setOption(
"max-degree",&maxdegree,
"Maximum degree of the GMRES polynomial.");
130 cmdp.
setOption(
"max-subspace",&maxsubspace,
"Maximum number of blocks the solver can use for the subspace.");
131 cmdp.
setOption(
"max-restarts",&maxrestarts,
"Maximum number of restarts allowed for GMRES solver.");
136 proc_verbose = verbose && (MyPID==0);
143 #ifndef HAVE_BELOS_TRIUTILS
144 std::cout <<
"This test requires Triutils. Please configure with --enable-triutils." << std::endl;
146 std::cout <<
"End Result: TEST FAILED" << std::endl;
157 info = readHB_newmat_double(filename.c_str(),&dim,&dim2,&nnz,
158 &colptr,&rowind,&dvals);
159 if (info == 0 || nnz < 0) {
161 std::cout <<
"Error reading '" << filename <<
"'" << std::endl;
162 std::cout <<
"End Result: TEST FAILED" << std::endl;
168 for (
int ii=0; ii<nnz; ii++) {
169 cvals[ii] = ST(dvals[ii*2],dvals[ii*2+1]);
181 MVT::MvRandom( *soln );
182 OPT::Apply( *A, *soln, *rhs );
183 MVT::MvInit( *soln, zero );
189 problem->setInitResVec( rhs );
190 bool set = problem->setProblem();
193 std::cout << std::endl <<
"ERROR: Belos::LinearProblem failed to set up correctly!" << std::endl;
201 maxiters = dim/blocksize - 1;
204 belosList.
set(
"Num Blocks", maxsubspace);
205 belosList.
set(
"Block Size", blocksize );
206 belosList.
set(
"Maximum Iterations", maxiters );
207 belosList.
set(
"Maximum Restarts", maxrestarts );
208 belosList.
set(
"Convergence Tolerance", tol );
213 belosList.
set(
"Output Frequency", frequency );
215 belosList.
set(
"Verbosity", verbosity );
218 polyList.
set(
"Maximum Degree", maxdegree );
219 polyList.
set(
"Polynomial Tolerance", polytol );
220 polyList.
set(
"Verbosity", verbosity );
221 polyList.
set(
"Random RHS", userandomrhs );
222 if ( outersolver !=
"" ) {
223 polyList.
set(
"Outer Solver", outersolver );
224 polyList.
set(
"Outer Solver Params", belosList );
245 std::cout << std::endl << std::endl;
246 std::cout <<
"Dimension of matrix: " << dim << std::endl;
247 std::cout <<
"Number of right-hand sides: " << numrhs << std::endl;
248 std::cout <<
"Block size used by solver: " << blocksize << std::endl;
249 std::cout <<
"Max number of Gmres iterations: " << maxiters << std::endl;
250 std::cout <<
"Relative residual tolerance: " << tol << std::endl;
251 std::cout << std::endl;
261 OPT::Apply( *A, *soln, *temp );
262 MVT::MvAddMv( one, *rhs, -one, *temp, *temp );
263 std::vector<MT> norm_num(numrhs), norm_denom(numrhs);
264 MVT::MvNorm( *temp, norm_num );
265 MVT::MvNorm( *rhs, norm_denom );
266 for (
int i=0; i<numrhs; ++i) {
268 std::cout <<
"Relative residual "<<i<<
" : " << norm_num[i] / norm_denom[i] << std::endl;
269 if ( norm_num[i] / norm_denom[i] > tol ) {
275 const std::vector<MT> residualLog = debugTest.
getLogResNorm();
276 if (numrhs==1 && proc_verbose && residualLog.size())
278 std::cout <<
"Absolute residual 2-norm [ " << residualLog.size() <<
" ] : ";
279 for (
unsigned int i=0; i<residualLog.size(); i++)
280 std::cout << residualLog[i] <<
" ";
281 std::cout << std::endl;
282 std::cout <<
"Final abs 2-norm / rhs 2-norm : " << residualLog[residualLog.size()-1] / norm_denom[0] << std::endl;
294 std::cout <<
"End Result: TEST PASSED" << std::endl;
297 std::cout <<
"End Result: TEST FAILED" << std::endl;
302 return ( success ? EXIT_SUCCESS : EXIT_FAILURE );
std::string Belos_Version()
int main(int argc, char *argv[])
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
const std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > & getLogResNorm() const
Returns the log of the absolute residual norm from the iteration.
A Belos::StatusTest debugging class for storing the absolute residual norms generated during a solve...
Belos::StatusTest debugging class for storing the absolute residual norms generated during a solve...
Traits class which defines basic operations on multivectors.
Simple example of a user's defined Belos::MultiVec class.
Alternative run-time polymorphic interface for operators.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void setOption(const char option_true[], const char option_false[], bool *option_val, const char documentation[]=NULL)
Declaration and definition of Belos::GmresPolySolMgr (hybrid block GMRES linear solver).
#define TEUCHOS_STANDARD_CATCH_STATEMENTS(VERBOSE, ERR_STREAM, SUCCESS_FLAG)
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
EParseCommandLineReturn parse(int argc, char *argv[], std::ostream *errout=&std::cerr) const
Hybrid block GMRES iterative linear solver.
ReturnType
Whether the Belos solve converged for all linear systems.
Interface for multivectors used by Belos' linear solvers.
Class which defines basic traits for the operator type.
Belos header file which uses auto-configuration information to include necessary C++ headers...
Simple example of a user's defined Belos::Operator class.