63 #ifdef HAVE_BELOS_TRIUTILS
64 #include "Trilinos_Util_iohb.h"
71 using namespace Teuchos;
73 int main(
int argc,
char *argv[]) {
75 typedef std::complex<double> ST;
77 typedef SCT::magnitudeType MT;
83 ST zero = SCT::zero();
96 bool norm_failure =
false;
97 bool proc_verbose =
false;
102 int maxrestarts = 15;
104 std::string
filename(
"mhd1280b.cua");
105 std::string ortho(
"ICGS");
109 cmdp.
setOption(
"verbose",
"quiet",&verbose,
"Print messages and results.");
110 cmdp.
setOption(
"debug",
"no-debug",&debug,
"Print debug messages.");
111 cmdp.
setOption(
"pseudo",
"regular",&pseudo,
"Use pseudo-block GMRES to solve the linear systems.");
112 cmdp.
setOption(
"frequency",&frequency,
"Solvers frequency for printing residuals (#iters).");
113 cmdp.
setOption(
"filename",&filename,
"Filename for Harwell-Boeing test matrix.");
114 cmdp.
setOption(
"tol",&tol,
"Relative residual tolerance used by GMRES solver.");
115 cmdp.
setOption(
"num-rhs",&numrhs,
"Number of right-hand sides to be solved for.");
116 cmdp.
setOption(
"num-restarts",&maxrestarts,
"Maximum number of restarts allowed for the GMRES solver.");
117 cmdp.
setOption(
"blocksize",&blocksize,
"Block size used by GMRES.");
118 cmdp.
setOption(
"subspace-length",&length,
"Maximum dimension of block-subspace used by GMRES solver.");
119 cmdp.
setOption(
"ortho-type",&ortho,
"Orthogonalization type, either DGKS, ICGS or IMGS");
124 proc_verbose = verbose && (MyPID==0);
131 #ifndef HAVE_BELOS_TRIUTILS
132 std::cout <<
"This test requires Triutils. Please configure with --enable-triutils." << std::endl;
134 std::cout <<
"End Result: TEST FAILED" << std::endl;
145 info = readHB_newmat_double(filename.c_str(),&dim,&dim2,&nnz,
146 &colptr,&rowind,&dvals);
147 if (info == 0 || nnz < 0) {
149 std::cout <<
"Error reading '" << filename <<
"'" << std::endl;
150 std::cout <<
"End Result: TEST FAILED" << std::endl;
156 for (
int ii=0; ii<nnz; ii++) {
157 cvals[ii] = ST(dvals[ii*2],dvals[ii*2+1]);
166 int maxits = dim/blocksize;
169 belosList.
set(
"Num Blocks", length );
170 belosList.
set(
"Block Size", blocksize );
171 belosList.
set(
"Maximum Iterations", maxits );
172 belosList.
set(
"Maximum Restarts", maxrestarts );
173 belosList.
set(
"Convergence Tolerance", tol );
174 belosList.
set(
"Orthogonalization", ortho );
180 belosList.
set(
"Verbosity", verbosity );
182 belosList.
set(
"Output Frequency", frequency );
193 MVT::MvRandom( *soln );
194 OPT::Apply( *A, *soln, *rhs );
195 MVT::MvInit( *soln, zero );
201 bool set = problem->setProblem();
204 std::cout << std::endl <<
"ERROR: Belos::LinearProblem failed to set up correctly!" << std::endl;
223 solver->setDebugStatusTest(
Teuchos::rcp(&debugTest,
false) );
229 std::cout << std::endl << std::endl;
230 std::cout <<
"Dimension of matrix: " << dim << std::endl;
231 std::cout <<
"Number of right-hand sides: " << numrhs << std::endl;
232 std::cout <<
"Block size used by solver: " << blocksize << std::endl;
233 std::cout <<
"Max number of Gmres iterations: " << maxits << std::endl;
234 std::cout <<
"Relative residual tolerance: " << tol << std::endl;
235 std::cout << std::endl;
245 OPT::Apply( *A, *soln, *temp );
246 MVT::MvAddMv( one, *rhs, -one, *temp, *temp );
247 std::vector<MT> norm_num(numrhs), norm_denom(numrhs);
248 MVT::MvNorm( *temp, norm_num );
249 MVT::MvNorm( *rhs, norm_denom );
250 for (
int i=0; i<numrhs; ++i) {
252 std::cout <<
"Relative residual "<<i<<
" : " << norm_num[i] / norm_denom[i] << std::endl;
253 if ( norm_num[i] / norm_denom[i] > tol ) {
259 const std::vector<MT> residualLog = debugTest.
getLogResNorm();
260 if (numrhs==1 && proc_verbose && residualLog.size())
262 std::cout <<
"Absolute residual 2-norm [ " << residualLog.size() <<
" ] : ";
263 for (
unsigned int i=0; i<residualLog.size(); i++)
264 std::cout << residualLog[i] <<
" ";
265 std::cout << std::endl;
266 std::cout <<
"Final abs 2-norm / rhs 2-norm : " << residualLog[residualLog.size()-1] / norm_denom[0] << std::endl;
278 std::cout <<
"End Result: TEST PASSED" << std::endl;
281 std::cout <<
"End Result: TEST FAILED" << std::endl;
286 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...
Interface to Block GMRES and Flexible GMRES.
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)
The Belos::BlockGmresSolMgr provides a solver manager for the BlockGmres linear solver.
void setOption(const char option_true[], const char option_false[], bool *option_val, const char documentation[]=NULL)
#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
Interface to standard and "pseudoblock" GMRES.
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.
The Belos::PseudoBlockGmresSolMgr provides a solver manager for the BlockGmres linear solver...
Belos header file which uses auto-configuration information to include necessary C++ headers...
Simple example of a user's defined Belos::Operator class.