50 #ifndef BELOS_TFQMR_ITER_HPP
51 #define BELOS_TFQMR_ITER_HPP
93 template <
class ScalarType,
class MV>
105 Rtilde(Teuchos::null),
D(Teuchos::null),
V(Teuchos::null)
126 template <
class ScalarType,
class MV,
class OP>
254 "Belos::TFQMRIter::setBlockSize(): Cannot use a block size that is not one.");
322 template <
class ScalarType,
class MV,
class OP>
338 stateStorageInitialized_(false),
345 template <
class ScalarType,
class MV,
class OP>
359 template <
class ScalarType,
class MV,
class OP>
362 if (!stateStorageInitialized_) {
368 stateStorageInitialized_ =
false;
379 "Belos::TFQMRIter::setStateSize(): linear problem does not specify multivectors to clone from.");
380 R_ = MVT::Clone( *tmp, 1 );
381 D_ = MVT::Clone( *tmp, 1 );
382 V_ = MVT::Clone( *tmp, 1 );
383 solnUpdate_ = MVT::Clone( *tmp, 1 );
387 stateStorageInitialized_ =
true;
394 template <
class ScalarType,
class MV,
class OP>
398 if (!stateStorageInitialized_)
402 "Belos::TFQMRIter::initialize(): Cannot initialize state storage!");
406 std::string errstr(
"Belos::TFQMRIter::initialize(): Specified multivectors must have a consistent length and width.");
414 std::invalid_argument, errstr );
416 std::invalid_argument, errstr );
419 if (newstate.
R != R_) {
421 MVT::Assign( *newstate.
R, *R_ );
427 W_ = MVT::CloneCopy( *R_ );
428 U_ = MVT::CloneCopy( *R_ );
429 Rtilde_ = MVT::CloneCopy( *R_ );
431 MVT::MvInit( *solnUpdate_ );
435 lp_->apply( *U_, *V_ );
436 AU_ = MVT::CloneCopy( *V_ );
441 MVT::MvNorm( *R_, tau_ );
442 MVT::MvDot( *R_, *Rtilde_, rho_old_ );
447 "Belos::TFQMRIter::initialize(): TFQMRIterState does not have initial residual.");
457 template <
class ScalarType,
class MV,
class OP>
463 if (initialized_ ==
false) {
472 ScalarType eta = STzero, beta = STzero;
481 "Belos::TFQMRIter::iterate(): current linear system has more than one vector!" );
487 while (stest_->checkStatus(
this) !=
Passed) {
489 for (
int iIter=0; iIter<2; iIter++)
497 MVT::MvDot( *V_, *Rtilde_, alpha_ );
498 alpha_[0] = rho_old_[0]/alpha_[0];
506 MVT::MvAddMv( STone, *W_, -alpha_[0], *AU_, *W_ );
513 MVT::MvAddMv( STone, *U_, (theta_[0]*theta_[0]/alpha_[0])*eta, *D_, *D_ );
524 MVT::MvAddMv( STone, *U_, -alpha_[0], *V_, *U_ );
527 lp_->apply( *U_, *AU_ );
534 MVT::MvNorm( *W_, theta_ );
535 theta_[0] /= tau_[0];
538 tau_[0] *= theta_[0]*cs_[0];
539 eta = cs_[0]*cs_[0]*alpha_[0];
546 MVT::MvAddMv( STone, *solnUpdate_, eta, *D_, *solnUpdate_ );
551 if ( tau_[0] == MTzero ) {
561 MVT::MvDot( *W_, *Rtilde_, rho_ );
562 beta = rho_[0]/rho_old_[0];
563 rho_old_[0] = rho_[0];
570 MVT::MvAddMv( STone, *W_, beta, *U_, *U_ );
573 MVT::MvAddMv( STone, *AU_, beta, *V_, *V_ );
576 lp_->apply( *U_, *AU_ );
579 MVT::MvAddMv( STone, *AU_, beta, *V_, *V_ );
592 #endif // BELOS_TFQMR_ITER_HPP
Teuchos::RCP< const MV > V
Collection of types and exceptions used within the Belos solvers.
This class implements the preconditioned transpose-free QMR algorithm for solving non-Hermitian linea...
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
Teuchos::ScalarTraits< ScalarType > SCT
SCT::magnitudeType MagnitudeType
OperatorTraits< ScalarType, MV, OP > OPT
TFQMRIterateFailure(const std::string &what_arg)
void setBlockSize(int blockSize)
Set the blocksize.
Class which manages the output and verbosity of the Belos solvers.
std::vector< ScalarType > alpha_
Teuchos::RCP< MV > Rtilde_
Pure virtual base class for defining the status testing capabilities of Belos.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Declaration of basic traits for the multivector type.
void initialize()
Initialize the solver with the initial vectors from the linear problem or random data.
TFQMRIter(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< OutputManager< ScalarType > > &printer, const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > &tester, Teuchos::ParameterList ¶ms)
Belos::TFQMRIter constructor.
virtual ~TFQMRIter()
Belos::TFQMRIter destructor.
Teuchos::RCP< const MV > R
The current residual basis.
Teuchos::RCP< const MV > getNativeResiduals(std::vector< MagnitudeType > *norms) const
Get the norms of the residuals native to the solver.
Teuchos::RCP< const MV > U
void resetNumIters(int iter=0)
Reset the iteration count.
Pure virtual base class which describes the basic interface to the linear solver iteration.
A pure virtual class for defining the status tests for the Belos iterative solvers.
Class which defines basic traits for the operator type.
Traits class which defines basic operations on multivectors.
Teuchos::RCP< MV > getCurrentUpdate() const
Get the current update to the linear system.
void initializeTFQMR(const TFQMRIterState< ScalarType, MV > &newstate)
Initialize the solver to an iterate, providing a complete state.
Teuchos::RCP< const MV > Rtilde
void setStateSize()
Method for initalizing the state storage needed by TFQMR.
void iterate()
This method performs TFQMR iterations until the status test indicates the need to stop or an error oc...
bool isInitialized()
States whether the solver has been initialized or not.
Teuchos::RCP< const MV > W
A linear system to solve, and its associated information.
Structure to contain pointers to TFQMRIter state variables.
Class which describes the linear problem to be solved by the iterative solver.
std::vector< MagnitudeType > cs_
bool stateStorageInitialized_
const Teuchos::RCP< StatusTest< ScalarType, MV, OP > > stest_
std::vector< MagnitudeType > tau_
int getNumIters() const
Get the current iteration count.
MultiVecTraits< ScalarType, MV > MVT
Teuchos::RCP< MV > solnUpdate_
std::vector< ScalarType > rho_old_
TFQMRIterateFailure is thrown when the TFQMRIter object is unable to compute the next iterate in the ...
int getBlockSize() const
Get the blocksize to be used by the iterative solver in solving this linear problem.
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > lp_
Teuchos::RCP< const MV > D
Belos header file which uses auto-configuration information to include necessary C++ headers...
std::vector< MagnitudeType > theta_
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Get a constant reference to the linear problem.
std::vector< ScalarType > rho_
const Teuchos::RCP< OutputManager< ScalarType > > om_
TFQMRIterState< ScalarType, MV > getState() const
Get the current state of the linear solver.