Belos Package Browser (Single Doxygen Collection)  Development
BelosRCGIter.hpp
Go to the documentation of this file.
1 //@HEADER
2 // ************************************************************************
3 //
4 // Belos: Block Linear Solvers Package
5 // Copyright 2004 Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 //@HEADER
41 
42 #ifndef BELOS_RCG_ITER_HPP
43 #define BELOS_RCG_ITER_HPP
44 
49 #include "BelosConfigDefs.hpp"
50 #include "BelosTypes.hpp"
51 
52 #include "BelosLinearProblem.hpp"
53 #include "BelosMatOrthoManager.hpp"
54 #include "BelosOutputManager.hpp"
55 #include "BelosStatusTest.hpp"
56 #include "BelosOperatorTraits.hpp"
57 #include "BelosMultiVecTraits.hpp"
58 
59 #include "Teuchos_LAPACK.hpp"
60 #include "Teuchos_SerialDenseMatrix.hpp"
61 #include "Teuchos_SerialDenseVector.hpp"
62 #include "Teuchos_ScalarTraits.hpp"
63 #include "Teuchos_ParameterList.hpp"
64 #include "Teuchos_TimeMonitor.hpp"
65 
66 // MLP Remove after debugging
67 #include <fstream>
68 #include <iomanip>
69 
81 namespace Belos {
82 
84 
85 
90  template <class ScalarType, class MV>
91  struct RCGIterState {
96  int curDim;
97 
99  Teuchos::RCP<MV> P;
100 
102  Teuchos::RCP<MV> Ap;
103 
105  Teuchos::RCP<MV> r;
106 
108  Teuchos::RCP<MV> z;
109 
111  bool existU;
112 
114  Teuchos::RCP<MV> U, AU;
115 
118  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > Alpha;
119  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > Beta;
120  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > D;
121  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > rTz_old;
122 
124  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > Delta;
125 
127  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > LUUTAU;
129  Teuchos::RCP<std::vector<int> > ipiv;
130 
131 
132  RCGIterState() : curDim(0), P(Teuchos::null), Ap(Teuchos::null), r(Teuchos::null),
133  z(Teuchos::null),
134  existU(false),
135  U(Teuchos::null), AU(Teuchos::null),
136  Alpha(Teuchos::null), Beta(Teuchos::null), D(Teuchos::null), rTz_old(Teuchos::null),
137  Delta(Teuchos::null), LUUTAU(Teuchos::null), ipiv(Teuchos::null)
138  {}
139  };
140 
142 
144 
145 
157  class RCGIterInitFailure : public BelosError {public:
158  RCGIterInitFailure(const std::string& what_arg) : BelosError(what_arg)
159  {}};
160 
167  class RCGIterFailure : public BelosError {public:
168  RCGIterFailure(const std::string& what_arg) : BelosError(what_arg)
169  {}};
170 
177  class RCGIterLAPACKFailure : public BelosError {public:
178  RCGIterLAPACKFailure(const std::string& what_arg) : BelosError(what_arg)
179  {}};
180 
182 
183 
184  template<class ScalarType, class MV, class OP>
185  class RCGIter : virtual public Iteration<ScalarType,MV,OP> {
186 
187  public:
188 
189  //
190  // Convenience typedefs
191  //
194  typedef Teuchos::ScalarTraits<ScalarType> SCT;
195  typedef typename SCT::magnitudeType MagnitudeType;
196 
198 
199 
207  RCGIter( const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > &problem,
208  const Teuchos::RCP<OutputManager<ScalarType> > &printer,
209  const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &tester,
210  Teuchos::ParameterList &params );
211 
213  virtual ~RCGIter() {};
215 
216 
218 
219 
232  void iterate();
233 
248  void initialize(RCGIterState<ScalarType,MV> &newstate);
249 
253  void initialize()
254  {
256  initialize(empty);
257  }
258 
260 
262 
263 
265  int getNumIters() const { return iter_; }
266 
268  void resetNumIters( int iter = 0 ) { iter_ = iter; }
269 
272  Teuchos::RCP<const MV> getNativeResiduals( std::vector<MagnitudeType> * /* norms */ ) const { return r_; }
273 
275 
280  Teuchos::RCP<MV> getCurrentUpdate() const { return Teuchos::null; }
281 
283  int getCurSubspaceDim() const {
284  if (!initialized_) return 0;
285  return curDim_;
286  };
287 
289  int getMaxSubspaceDim() const { return numBlocks_+1; }
290 
292 
293 
295 
296 
298  const LinearProblem<ScalarType,MV,OP>& getProblem() const { return *lp_; }
299 
301  int getNumBlocks() const { return numBlocks_; }
302 
304  void setNumBlocks(int numBlocks) { setSize( recycleBlocks_, numBlocks ); };
305 
307  int getRecycledBlocks() const { return recycleBlocks_; }
308 
310  void setRecycledBlocks(int recycleBlocks) { setSize( recycleBlocks, numBlocks_ ); };
311 
313  int getBlockSize() const { return 1; }
314 
316  void setBlockSize(int blockSize) {
317  TEUCHOS_TEST_FOR_EXCEPTION(blockSize!=1,std::invalid_argument,
318  "Belos::RCGIter::setBlockSize(): Cannot use a block size that is not one.");
319  }
320 
322  void setSize( int recycleBlocks, int numBlocks );
323 
325  bool isInitialized() { return initialized_; }
326 
328 
329  private:
330 
331  //
332  // Internal methods
333  //
334 
335  //
336  // Classes input through constructor that define the linear problem to be solved.
337  //
338  const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > lp_;
339  const Teuchos::RCP<OutputManager<ScalarType> > om_;
340  const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > stest_;
341 
342  //
343  // Algorithmic parameters
344  //
345  // numBlocks_ is the size of the allocated space for the Krylov basis, in blocks.
347 
348  // recycleBlocks_ is the size of the allocated space for the recycled subspace, in blocks.
350 
351  //
352  // Current solver state
353  //
354  // initialized_ specifies that the basis vectors have been initialized and the iterate() routine
355  // is capable of running; _initialize is controlled by the initialize() member method
356  // For the implications of the state of initialized_, please see documentation for initialize()
358 
359  // Current subspace dimension, and number of iterations performed.
361 
362  //
363  // State Storage
364  //
365  // Search vectors
366  Teuchos::RCP<MV> P_;
367  //
368  // A times current search vector
369  Teuchos::RCP<MV> Ap_;
370  //
371  // Residual vector
372  Teuchos::RCP<MV> r_;
373  //
374  // Preconditioned residual
375  Teuchos::RCP<MV> z_;
376  //
377  // Flag to indicate that the recycle space should be used
378  bool existU_;
379  // Recycled subspace and its image
380  Teuchos::RCP<MV> U_, AU_;
381  //
382  // Coefficients arising in RCG iteration
383  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > Alpha_,Beta_,D_;
384  //
385  // Solutions to local least-squares problems
386  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > Delta_;
387  //
388  // The LU factorization of the matrix U^T A U
389  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > LUUTAU_;
390  //
391  // Data from LU factorization of UTAU
392  Teuchos::RCP<std::vector<int> > ipiv_;
393  //
394  // The scalar r'*z
395  Teuchos::RCP<Teuchos::SerialDenseMatrix<int,ScalarType> > rTz_old_;
396  };
397 
399  // Constructor.
400  template<class ScalarType, class MV, class OP>
402  const Teuchos::RCP<OutputManager<ScalarType> > &printer,
403  const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > &tester,
404  Teuchos::ParameterList &params ):
405  lp_(problem),
406  om_(printer),
407  stest_(tester),
408  numBlocks_(0),
409  recycleBlocks_(0),
410  initialized_(false),
411  curDim_(0),
412  iter_(0),
413  existU_(false)
414  {
415  // Get the maximum number of blocks allowed for this Krylov subspace
416  TEUCHOS_TEST_FOR_EXCEPTION(!params.isParameter("Num Blocks"), std::invalid_argument,
417  "Belos::RCGIter::constructor: mandatory parameter \"Num Blocks\" is not specified.");
418  int nb = Teuchos::getParameter<int>(params, "Num Blocks");
419 
420  TEUCHOS_TEST_FOR_EXCEPTION(!params.isParameter("Recycled Blocks"), std::invalid_argument,
421  "Belos::RCGIter::constructor: mandatory parameter \"Recycled Blocks\" is not specified.");
422  int rb = Teuchos::getParameter<int>(params, "Recycled Blocks");
423 
424  // Set the number of blocks and allocate data
425  setSize( rb, nb );
426  }
427 
429  // Set the block size and make necessary adjustments.
430  template <class ScalarType, class MV, class OP>
431  void RCGIter<ScalarType,MV,OP>::setSize( int recycleBlocks, int numBlocks )
432  {
433 
434  TEUCHOS_TEST_FOR_EXCEPTION(numBlocks <= 0, std::invalid_argument, "Belos::RCGIter::setSize() was passed a non-positive argument for \"Num Blocks\".");
435  TEUCHOS_TEST_FOR_EXCEPTION(recycleBlocks <= 0, std::invalid_argument, "Belos::RCGIter::setSize() was passed a non-positive argument for \"Recycled Blocks\".");
436  TEUCHOS_TEST_FOR_EXCEPTION(recycleBlocks >= numBlocks, std::invalid_argument, "Belos::RCGIter::setSize() the number of recycled blocks is larger than the allowable subspace.");
437 
438  numBlocks_ = numBlocks;
439  recycleBlocks_ = recycleBlocks;
440 
441  }
442 
444  // Initialize this iteration object
445  template <class ScalarType, class MV, class OP>
447  {
448 
449  if (newstate.P != Teuchos::null &&
450  newstate.Ap != Teuchos::null &&
451  newstate.r != Teuchos::null &&
452  newstate.z != Teuchos::null &&
453  newstate.U != Teuchos::null &&
454  newstate.AU != Teuchos::null &&
455  newstate.Alpha != Teuchos::null &&
456  newstate.Beta != Teuchos::null &&
457  newstate.D != Teuchos::null &&
458  newstate.Delta != Teuchos::null &&
459  newstate.LUUTAU != Teuchos::null &&
460  newstate.ipiv != Teuchos::null &&
461  newstate.rTz_old != Teuchos::null) {
462 
463  curDim_ = newstate.curDim;
464  P_ = newstate.P;
465  Ap_ = newstate.Ap;
466  r_ = newstate.r;
467  z_ = newstate.z;
468  existU_ = newstate.existU;
469  U_ = newstate.U;
470  AU_ = newstate.AU;
471  Alpha_ = newstate.Alpha;
472  Beta_ = newstate.Beta;
473  D_ = newstate.D;
474  Delta_ = newstate.Delta;
475  LUUTAU_ = newstate.LUUTAU;
476  ipiv_ = newstate.ipiv;
477  rTz_old_ = newstate.rTz_old;
478  }
479  else {
480 
481  TEUCHOS_TEST_FOR_EXCEPTION(newstate.P == Teuchos::null,std::invalid_argument,
482  "Belos::RCGIter::initialize(): RCGIterState does not have P initialized.");
483 
484  TEUCHOS_TEST_FOR_EXCEPTION(newstate.Ap == Teuchos::null,std::invalid_argument,
485  "Belos::RCGIter::initialize(): RCGIterState does not have Ap initialized.");
486 
487  TEUCHOS_TEST_FOR_EXCEPTION(newstate.r == Teuchos::null,std::invalid_argument,
488  "Belos::RCGIter::initialize(): RCGIterState does not have r initialized.");
489 
490  TEUCHOS_TEST_FOR_EXCEPTION(newstate.z == Teuchos::null,std::invalid_argument,
491  "Belos::RCGIter::initialize(): RCGIterState does not have z initialized.");
492 
493  TEUCHOS_TEST_FOR_EXCEPTION(newstate.U == Teuchos::null,std::invalid_argument,
494  "Belos::RCGIter::initialize(): RCGIterState does not have U initialized.");
495 
496  TEUCHOS_TEST_FOR_EXCEPTION(newstate.AU == Teuchos::null,std::invalid_argument,
497  "Belos::RCGIter::initialize(): RCGIterState does not have AU initialized.");
498 
499  TEUCHOS_TEST_FOR_EXCEPTION(newstate.Alpha == Teuchos::null,std::invalid_argument,
500  "Belos::RCGIter::initialize(): RCGIterState does not have Alpha initialized.");
501 
502  TEUCHOS_TEST_FOR_EXCEPTION(newstate.Beta == Teuchos::null,std::invalid_argument,
503  "Belos::RCGIter::initialize(): RCGIterState does not have Beta initialized.");
504 
505  TEUCHOS_TEST_FOR_EXCEPTION(newstate.D == Teuchos::null,std::invalid_argument,
506  "Belos::RCGIter::initialize(): RCGIterState does not have D initialized.");
507 
508  TEUCHOS_TEST_FOR_EXCEPTION(newstate.Delta == Teuchos::null,std::invalid_argument,
509  "Belos::RCGIter::initialize(): RCGIterState does not have Delta initialized.");
510 
511  TEUCHOS_TEST_FOR_EXCEPTION(newstate.LUUTAU == Teuchos::null,std::invalid_argument,
512  "Belos::RCGIter::initialize(): RCGIterState does not have LUUTAU initialized.");
513 
514  TEUCHOS_TEST_FOR_EXCEPTION(newstate.ipiv == Teuchos::null,std::invalid_argument,
515  "Belos::RCGIter::initialize(): RCGIterState does not have ipiv initialized.");
516 
517  TEUCHOS_TEST_FOR_EXCEPTION(newstate.rTz_old == Teuchos::null,std::invalid_argument,
518  "Belos::RCGIter::initialize(): RCGIterState does not have rTz_old initialized.");
519 
520  }
521 
522  // the solver is initialized
523  initialized_ = true;
524 
525  }
526 
528  // Iterate until the status test informs us we should stop.
529  template <class ScalarType, class MV, class OP>
531  {
532  TEUCHOS_TEST_FOR_EXCEPTION( initialized_ == false, RCGIterFailure,
533  "Belos::RCGIter::iterate(): RCGIter class not initialized." );
534 
535  // We'll need LAPACK
536  Teuchos::LAPACK<int,ScalarType> lapack;
537 
538  // Create convenience variables for zero and one.
539  ScalarType one = Teuchos::ScalarTraits<ScalarType>::one();
540  ScalarType zero = Teuchos::ScalarTraits<ScalarType>::zero();
541 
542  // Allocate memory for scalars
543  std::vector<int> index(1);
544  Teuchos::SerialDenseMatrix<int,ScalarType> pAp(1,1), rTz(1,1);
545 
546  // Get the current solution std::vector.
547  Teuchos::RCP<MV> cur_soln_vec = lp_->getCurrLHSVec();
548 
549  // Check that the current solution std::vector only has one column.
550  TEUCHOS_TEST_FOR_EXCEPTION( MVT::GetNumberVecs(*cur_soln_vec) != 1, RCGIterFailure,
551  "Belos::RCGIter::iterate(): current linear system has more than one std::vector!" );
552 
553  // Compute the current search dimension.
554  int searchDim = numBlocks_+1;
555 
556  // index of iteration within current cycle
557  int i_ = 0;
558 
560  // iterate until the status test tells us to stop.
561  //
562  // also break if our basis is full
563  //
564  Teuchos::RCP<const MV> p_ = Teuchos::null;
565  Teuchos::RCP<MV> pnext_ = Teuchos::null;
566  while (stest_->checkStatus(this) != Passed && curDim_+1 <= searchDim) {
567 
568  // Ap = A*p;
569  index.resize( 1 );
570  index[0] = i_;
571  p_ = MVT::CloneView( *P_, index );
572  lp_->applyOp( *p_, *Ap_ );
573 
574  // d = p'*Ap;
575  MVT::MvTransMv( one, *p_, *Ap_, pAp );
576  (*D_)(i_,0) = pAp(0,0);
577 
578  // alpha = rTz_old / pAp
579  (*Alpha_)(i_,0) = (*rTz_old_)(0,0) / pAp(0,0);
580 
581  // Check that alpha is a positive number
582  TEUCHOS_TEST_FOR_EXCEPTION( SCT::real(pAp(0,0)) <= zero, RCGIterFailure, "Belos::RCGIter::iterate(): non-positive value for p^H*A*p encountered!" );
583 
584  // x = x + (alpha * p);
585  MVT::MvAddMv( one, *cur_soln_vec, (*Alpha_)(i_,0), *p_, *cur_soln_vec );
586  lp_->updateSolution();
587 
588  // r = r - (alpha * Ap);
589  MVT::MvAddMv( one, *r_, -(*Alpha_)(i_,0), *Ap_, *r_ );
590 
591  std::vector<MagnitudeType> norm(1);
592  MVT::MvNorm( *r_, norm );
593 //printf("i = %i\tnorm(r) = %e\n",i_,norm[0]);
594 
595  // z = M\r
596  if ( lp_->getLeftPrec() != Teuchos::null ) {
597  lp_->applyLeftPrec( *r_, *z_ );
598  }
599  else if ( lp_->getRightPrec() != Teuchos::null ) {
600  lp_->applyRightPrec( *r_, *z_ );
601  }
602  else {
603  z_ = r_;
604  }
605 
606  // rTz_new = r'*z;
607  MVT::MvTransMv( one, *r_, *z_, rTz );
608 
609  // beta = rTz_new/rTz_old;
610  (*Beta_)(i_,0) = rTz(0,0) / (*rTz_old_)(0,0);
611 
612  // rTz_old = rTz_new;
613  (*rTz_old_)(0,0) = rTz(0,0);
614 
615  // get pointer for next p
616  index.resize( 1 );
617  index[0] = i_+1;
618  pnext_ = MVT::CloneViewNonConst( *P_, index );
619 
620  if (existU_) {
621  // mu = UTAU \ (AU'*z);
622  Teuchos::SerialDenseMatrix<int,ScalarType> mu( Teuchos::View, *Delta_, recycleBlocks_, 1, 0, i_ );
623  MVT::MvTransMv( one, *AU_, *z_, mu );
624  char TRANS = 'N';
625  int info;
626  lapack.GETRS( TRANS, recycleBlocks_, 1, LUUTAU_->values(), LUUTAU_->stride(), &(*ipiv_)[0], mu.values(), mu.stride(), &info );
627  TEUCHOS_TEST_FOR_EXCEPTION(info != 0, RCGIterLAPACKFailure,
628  "Belos::RCGIter::solve(): LAPACK GETRS failed to compute a solution.");
629  // p = -(U*mu) + (beta*p) + z (in two steps)
630  // p = (beta*p) + z;
631  MVT::MvAddMv( (*Beta_)(i_,0), *p_, one, *z_, *pnext_ );
632  // pnext = -(U*mu) + (one)*pnext;
633  MVT::MvTimesMatAddMv( -one, *U_, mu, one, *pnext_ );
634  }
635  else {
636  // p = (beta*p) + z;
637  MVT::MvAddMv( (*Beta_)(i_,0), *p_, one, *z_, *pnext_ );
638  }
639 
640  // Done with this view; release pointer
641  p_ = Teuchos::null;
642  pnext_ = Teuchos::null;
643 
644  // increment iteration count and dimension index
645  i_++;
646  iter_++;
647  curDim_++;
648 
649  } // end while (statusTest == false)
650 
651  }
652 
653 } // end Belos namespace
654 
655 #endif /* BELOS_RCG_ITER_HPP */
const Teuchos::RCP< OutputManager< ScalarType > > om_
Teuchos::RCP< MV > r
The current residual.
Collection of types and exceptions used within the Belos solvers.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > rTz_old
Belos&#39;s basic output manager for sending information of select verbosity levels to the appropriate ou...
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Alpha_
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > rTz_old_
Teuchos::RCP< MV > P
The current Krylov basis.
Teuchos::RCP< std::vector< int > > ipiv
Data from LU factorization of U^T A U.
Class which manages the output and verbosity of the Belos solvers.
int getNumBlocks() const
Get the maximum number of blocks used by the iterative solver in solving this linear problem...
void setNumBlocks(int numBlocks)
Set the maximum number of blocks used by the iterative solver.
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Get a constant reference to the linear problem.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Beta_
Teuchos::RCP< MV > AU
RCGIterFailure is thrown when the RCGIter object is unable to compute the next iterate in the RCGIter...
int getBlockSize() const
Get the blocksize to be used by the iterative solver in solving this linear problem.
int getMaxSubspaceDim() const
Get the maximum dimension allocated for the search subspace.
Pure virtual base class for defining the status testing capabilities of Belos.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > D_
Teuchos::RCP< MV > U
The recycled subspace and its image.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > LUUTAU
The LU factorization of the matrix U^T A U.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Beta
Declaration of basic traits for the multivector type.
This class implements the RCG iteration, where a single-std::vector Krylov subspace is constructed...
const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > lp_
OperatorTraits< ScalarType, MV, OP > OPT
A pure virtual class for defining the status tests for the Belos iterative solvers.
Class which defines basic traits for the operator type.
MultiVecTraits< ScalarType, MV > MVT
virtual ~RCGIter()
Destructor.
RCGIterLAPACKFailure(const std::string &what_arg)
int getNumIters() const
Get the current iteration count.
Traits class which defines basic operations on multivectors.
void setSize(int recycleBlocks, int numBlocks)
Set the maximum number of blocks used by the iterative solver and the number of recycled vectors...
const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > stest_
Structure to contain pointers to RCGIter state variables.
bool existU
Flag to indicate the recycle space should be used.
Teuchos::RCP< MV > z
The current preconditioned residual.
Teuchos::RCP< MV > Ap
A times current search vector.
RCGIterInitFailure(const std::string &what_arg)
RCGIter(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< OutputManager< ScalarType > > &printer, const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, Teuchos::ParameterList &params)
RCGIter constructor with linear problem, solver utilities, and parameter list of solver options...
Teuchos::RCP< MV > U_
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
RCGIterFailure(const std::string &what_arg)
Teuchos::RCP< MV > P_
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Alpha
Coefficients arising in RCG iteration.
Teuchos::RCP< MV > getCurrentUpdate() const
Get the current update to the linear system.
int getRecycledBlocks() const
Get the maximum number of recycled blocks used by the iterative solver in solving this linear problem...
Teuchos::RCP< MV > Ap_
RCGIterInitFailure is thrown when the RCGIter object is unable to generate an initial iterate in the ...
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > LUUTAU_
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Delta_
int curDim
The current dimension of the reduction.
void setRecycledBlocks(int recycleBlocks)
Set the maximum number of recycled blocks used by the iterative solver.
void initialize()
Initialize the solver with the initial vectors from the linear problem or random data.
Teuchos::RCP< MV > AU_
int getCurSubspaceDim() const
Get the dimension of the search subspace used to generate the current solution to the linear problem...
Class which defines basic traits for the operator type.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > D
void iterate()
This method performs RCG iterations until the status test indicates the need to stop or an error occu...
Parent class to all Belos exceptions.
Definition: BelosTypes.hpp:60
void resetNumIters(int iter=0)
Reset the iteration count.
Teuchos::RCP< std::vector< int > > ipiv_
RCGIterLAPACKFailure is thrown when a nonzero return value is passed back from an LAPACK routine...
void setBlockSize(int blockSize)
Set the blocksize.
Belos header file which uses auto-configuration information to include necessary C++ headers...
Teuchos::ScalarTraits< ScalarType > SCT
Teuchos::RCP< MV > r_
bool isInitialized()
States whether the solver has been initialized or not.
Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > Delta
Solutions to local least-squares problems.
Teuchos::RCP< MV > z_
Teuchos::RCP< const MV > getNativeResiduals(std::vector< MagnitudeType > *) const
Get the norms of the residuals native to the solver.
Templated virtual class for providing orthogonalization/orthonormalization methods with matrix-based ...
SCT::magnitudeType MagnitudeType