87 typedef MultiVecTraits<ScalarType,MV> MVT;
88 typedef OperatorTraits<ScalarType,MV,OP> OPT;
89 typedef Teuchos::ScalarTraits<ScalarType> SCT;
90 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
91 typedef Teuchos::ScalarTraits<MagnitudeType> MT;
116 RTRSolMgr(
const Teuchos::RCP<Eigenproblem<ScalarType,MV,OP> > &problem,
117 Teuchos::ParameterList &pl );
136 Teuchos::Array<Teuchos::RCP<Teuchos::Time> >
getTimers()
const {
137 return Teuchos::tuple(_timerSolve);
163 Teuchos::RCP<Eigenproblem<ScalarType,MV,OP> > problem_;
168 MagnitudeType convtol_;
176 Teuchos::RCP<Teuchos::Time> _timerSolve;
177 Teuchos::RCP<OutputManager<ScalarType> > printer_;
178 Teuchos::ParameterList& pl_;
185 const Teuchos::RCP<Eigenproblem<ScalarType,MV,OP> > &problem,
186 Teuchos::ParameterList &pl ) :
190 convtol_(MT::prec()),
194#ifdef ANASAZI_TEUCHOS_TIME_MONITOR
195 _timerSolve(Teuchos::TimeMonitor::getNewTimer(
"Anasazi: RTRSolMgr::solve()")),
199 TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument,
"Problem not given to solver manager.");
200 TEUCHOS_TEST_FOR_EXCEPTION(!problem_->isProblemSet(), std::invalid_argument,
"Problem not set.");
201 TEUCHOS_TEST_FOR_EXCEPTION(!problem_->isHermitian(), std::invalid_argument,
"Problem not symmetric.");
202 TEUCHOS_TEST_FOR_EXCEPTION(problem_->getInitVec() == Teuchos::null,std::invalid_argument,
"Problem does not contain initial vectors to clone from.");
206 whch_ = pl_.get(
"Which",
"SR");
207 TEUCHOS_TEST_FOR_EXCEPTION(whch_ !=
"SR" && whch_ !=
"LR",
208 std::invalid_argument,
"Anasazi::RTRSolMgr: Invalid sorting string. RTR solvers compute only LR or SR.");
211 convtol_ = pl_.get(
"Convergence Tolerance",convtol_);
212 relconvtol_ = pl_.get(
"Relative Convergence Tolerance",relconvtol_);
213 strtmp = pl_.get(
"Convergence Norm",std::string(
"2"));
215 convNorm_ = RES_2NORM;
217 else if (strtmp ==
"M") {
218 convNorm_ = RES_ORTH;
221 TEUCHOS_TEST_FOR_EXCEPTION(
true, std::invalid_argument,
222 "Anasazi::RTRSolMgr: Invalid Convergence Norm.");
227 maxIters_ = pl_.get(
"Maximum Iterations",maxIters_);
230 skinny_ = pl_.get(
"Skinny Solver",skinny_);
233 numICGS_ = pl_.get(
"Num ICGS",2);
236 blkSize_ = pl_.get(
"Block Size", problem_->getNEV());
240 int osProc = pl.get(
"Output Processor", 0);
243 Teuchos::RCP<Teuchos::FancyOStream> osp;
245 if (pl.isParameter(
"Output Stream")) {
246 osp = Teuchos::getParameter<Teuchos::RCP<Teuchos::FancyOStream> >(pl,
"Output Stream");
253 if (pl_.isParameter(
"Verbosity")) {
254 if (Teuchos::isParameterType<int>(pl_,
"Verbosity")) {
255 verbosity = pl_.get(
"Verbosity", verbosity);
257 verbosity = (int)Teuchos::getParameter<Anasazi::MsgType>(pl_,
"Verbosity");
260 printer_ = Teuchos::rcp(
new OutputManager<ScalarType>(verbosity,osp) );
274 const int nev = problem_->getNEV();
280 Teuchos::RCP<Teuchos::FancyOStream>
281 out = Teuchos::getFancyOStream(Teuchos::rcpFromRef(printer_->stream(
Debug)));
282 out->setShowAllFrontMatter(
false).setShowProcRank(
true);
283 *out <<
"Entering Anasazi::RTRSolMgr::solve()\n";
288 Teuchos::RCP<BasicSort<MagnitudeType> > sorter = Teuchos::rcp(
new BasicSort<MagnitudeType>(whch_) );
293 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > maxtest;
294 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > ordertest;
295 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > combotest;
296 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convtest;
299 maxtest = Teuchos::rcp(
new StatusTestMaxIters<ScalarType,MV,OP>(maxIters_) );
302 maxtest = Teuchos::null;
306 convtest = Teuchos::rcp(
new StatusTestResNorm<ScalarType,MV,OP>(convtol_,nev,convNorm_,relconvtol_) );
307 ordertest = Teuchos::rcp(
new StatusTestWithOrdering<ScalarType,MV,OP>(convtest,sorter,nev) );
310 Teuchos::Array<Teuchos::RCP<StatusTest<ScalarType,MV,OP> > > alltests;
311 alltests.push_back(ordertest);
312 if (maxtest != Teuchos::null) alltests.push_back(maxtest);
316 Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputtest;
317 if ( printer_->isVerbosity(
Debug) ) {
318 outputtest = Teuchos::rcp(
new StatusTestOutput<ScalarType,MV,OP>( printer_,combotest,1,
Passed+
Failed+
Undefined ) );
321 outputtest = Teuchos::rcp(
new StatusTestOutput<ScalarType,MV,OP>( printer_,combotest,1,
Passed ) );
326 Teuchos::RCP<ICGSOrthoManager<ScalarType,MV,OP> > ortho
327 = Teuchos::rcp(
new ICGSOrthoManager<ScalarType,MV,OP>(problem_->getM(),numICGS_) );
332 bool leftMost =
true;
333 if (whch_ ==
"LR" || whch_ ==
"LM") {
336 pl_.set<
bool>(
"Leftmost",leftMost);
337 Teuchos::RCP<RTRBase<ScalarType,MV,OP> > rtr_solver;
338 if (skinny_ ==
false) {
340 rtr_solver = Teuchos::rcp(
new IRTR<ScalarType,MV,OP>(problem_,sorter,printer_,outputtest,ortho,pl_) );
344 rtr_solver = Teuchos::rcp(
new SIRTR<ScalarType,MV,OP>(problem_,sorter,printer_,outputtest,ortho,pl_) );
347 Teuchos::RCP< const MV > probauxvecs = problem_->getAuxVecs();
348 if (probauxvecs != Teuchos::null) {
349 rtr_solver->setAuxVecs( Teuchos::tuple< Teuchos::RCP<const MV> >(probauxvecs) );
352 TEUCHOS_TEST_FOR_EXCEPTION(rtr_solver->getBlockSize() < problem_->getNEV(),std::logic_error,
353 "Anasazi::RTRSolMgr requires block size >= requested number of eigenvalues.");
356 Teuchos::RCP<MV> foundvecs;
357 std::vector<MagnitudeType> foundvals;
361#ifdef ANASAZI_TEUCHOS_TIME_MONITOR
362 Teuchos::TimeMonitor slvtimer(*_timerSolve);
364 rtr_solver->iterate();
365 numIters_ = rtr_solver->getNumIters();
367 catch (
const std::exception &e) {
369 printer_->stream(
Anasazi::Errors) <<
"Exception: " << e.what() << endl;
370 Eigensolution<ScalarType,MV> sol;
372 problem_->setSolution(sol);
377 if (convtest->getStatus() ==
Passed || (maxtest != Teuchos::null && maxtest->getStatus() ==
Passed))
379 int num = convtest->howMany();
381 std::vector<int> ind = convtest->whichVecs();
383 foundvecs = MVT::CloneCopy(*rtr_solver->getRitzVectors(),ind);
385 foundvals.resize(num);
386 std::vector<Value<ScalarType> > all = rtr_solver->getRitzValues();
387 for (
int i=0; i<num; i++) {
388 foundvals[i] = all[ind[i]].realpart;
394 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"Anasazi::RTRSolMgr::solve(): solver returned without satisfy status test.");
398 Eigensolution<ScalarType,MV> sol;
399 sol.numVecs = numfound;
400 sol.Evecs = foundvecs;
401 sol.Espace = sol.Evecs;
402 sol.Evals.resize(sol.numVecs);
403 for (
int i=0; i<sol.numVecs; i++) {
404 sol.Evals[i].realpart = foundvals[i];
407 sol.index.resize(numfound,0);
410 rtr_solver->currentStatus(printer_->stream(
FinalSummary));
413#ifdef ANASAZI_TEUCHOS_TIME_MONITOR
415 Teuchos::TimeMonitor::summarize( printer_->stream(
TimingDetails ) );
420 problem_->setSolution(sol);
421 printer_->stream(
Debug) <<
"Returning " << sol.numVecs <<
" eigenpairs to eigenproblem." << endl;