50 #ifndef BELOS_TFQMR_ITER_HPP 51 #define BELOS_TFQMR_ITER_HPP 71 #include "Teuchos_BLAS.hpp" 72 #include "Teuchos_SerialDenseMatrix.hpp" 73 #include "Teuchos_SerialDenseVector.hpp" 74 #include "Teuchos_ScalarTraits.hpp" 75 #include "Teuchos_ParameterList.hpp" 76 #include "Teuchos_TimeMonitor.hpp" 95 template <
class ScalarType,
class MV>
99 Teuchos::RCP<const MV>
R;
100 Teuchos::RCP<const MV>
W;
101 Teuchos::RCP<const MV>
U;
103 Teuchos::RCP<const MV>
D;
104 Teuchos::RCP<const MV>
V;
107 Rtilde(Teuchos::null), D(Teuchos::null), V(Teuchos::null)
143 template <
class ScalarType,
class MV,
class OP>
151 typedef Teuchos::ScalarTraits<ScalarType>
SCT;
161 Teuchos::ParameterList ¶ms );
223 initializeTFQMR(empty);
241 state.solnUpdate = solnUpdate_;
259 Teuchos::RCP<const MV> getNativeResiduals( std::vector<MagnitudeType> *norms )
const;
281 TEUCHOS_TEST_FOR_EXCEPTION(blockSize!=1,std::invalid_argument,
282 "Belos::TFQMRIter::setBlockSize(): Cannot use a block size that is not one.");
302 const Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > lp_;
303 const Teuchos::RCP<OutputManager<ScalarType> > om_;
304 const Teuchos::RCP<StatusTest<ScalarType,MV,OP> > stest_;
311 Teuchos::SerialDenseMatrix<int,ScalarType> alpha_, rho_, rho_old_;
312 std::vector<MagnitudeType> tau_, cs_, theta_;
325 bool stateStorageInitialized_;
335 Teuchos::RCP<MV> U_, AU_;
336 Teuchos::RCP<MV> Rtilde_;
339 Teuchos::RCP<MV> solnUpdate_;
349 template <
class ScalarType,
class MV,
class OP>
353 Teuchos::ParameterList ¶ms
365 stateStorageInitialized_(false),
372 template <
class ScalarType,
class MV,
class OP>
373 Teuchos::RCP<const MV>
376 MagnitudeType one = Teuchos::ScalarTraits<MagnitudeType>::one();
378 (*normvec)[0] = Teuchos::ScalarTraits<MagnitudeType>::squareroot( iter_ + one )*tau_[0];
380 return Teuchos::null;
386 template <
class ScalarType,
class MV,
class OP>
389 if (!stateStorageInitialized_) {
392 Teuchos::RCP<const MV> lhsMV = lp_->getLHS();
393 Teuchos::RCP<const MV> rhsMV = lp_->getRHS();
394 if (lhsMV == Teuchos::null && rhsMV == Teuchos::null) {
395 stateStorageInitialized_ =
false;
402 if (R_ == Teuchos::null) {
404 Teuchos::RCP<const MV> tmp = ( (rhsMV!=Teuchos::null)? rhsMV: lhsMV );
405 TEUCHOS_TEST_FOR_EXCEPTION(tmp == Teuchos::null,std::invalid_argument,
406 "Belos::TFQMRIter::setStateSize(): linear problem does not specify multivectors to clone from.");
415 stateStorageInitialized_ =
true;
422 template <
class ScalarType,
class MV,
class OP>
426 if (!stateStorageInitialized_)
429 TEUCHOS_TEST_FOR_EXCEPTION(!stateStorageInitialized_,std::invalid_argument,
430 "Belos::TFQMRIter::initialize(): Cannot initialize state storage!");
434 std::string errstr(
"Belos::TFQMRIter::initialize(): Specified multivectors must have a consistent length and width.");
437 const ScalarType one = Teuchos::ScalarTraits<ScalarType>::one();
438 const ScalarType STzero = Teuchos::ScalarTraits<ScalarType>::zero();
439 const MagnitudeType MTzero = Teuchos::ScalarTraits<MagnitudeType>::zero();
441 if (newstate.
R != Teuchos::null) {
444 std::invalid_argument, errstr );
446 std::invalid_argument, errstr );
449 if (newstate.
R != R_) {
465 lp_->apply( *U_, *V_ );
476 TEUCHOS_TEST_FOR_EXCEPTION(newstate.
R == Teuchos::null,std::invalid_argument,
477 "Belos::TFQMRIter::initialize(): TFQMRIterState does not have initial residual.");
487 template <
class ScalarType,
class MV,
class OP>
493 if (initialized_ ==
false) {
498 const ScalarType STone = Teuchos::ScalarTraits<ScalarType>::one();
499 const MagnitudeType MTone = Teuchos::ScalarTraits<MagnitudeType>::one();
500 const ScalarType STzero = Teuchos::ScalarTraits<ScalarType>::zero();
501 ScalarType eta = STzero, beta = STzero;
506 Teuchos::RCP<MV> cur_soln_vec = lp_->getCurrLHSVec();
510 "Belos::TFQMRIter::iterate(): current linear system has more than one vector!" );
516 while (stest_->checkStatus(
this) !=
Passed) {
525 alpha_(0,0) = rho_old_(0,0)/alpha_(0,0);
540 MVT::MvAddMv( STone, *U_, (theta_[0]*theta_[0]/alpha_(0,0))*eta, *D_, *D_ );
554 lp_->apply( *U_, *AU_ );
562 theta_[0] /= tau_[0];
564 cs_[0] = MTone / Teuchos::ScalarTraits<MagnitudeType>::squareroot(MTone + theta_[0]*theta_[0]);
565 tau_[0] *= theta_[0]*cs_[0];
566 eta = cs_[0]*cs_[0]*alpha_(0,0);
574 MVT::MvAddMv( STone, *solnUpdate_, eta, *D_, *solnUpdate_ );
583 beta = rho_(0,0)/rho_old_(0,0);
584 rho_old_(0,0) = rho_(0,0);
597 lp_->apply( *U_, *AU_ );
611 #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.
Pure virtual base class for defining the status testing capabilities of Belos.
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.
static void MvInit(MV &mv, const ScalarType alpha=Teuchos::ScalarTraits< ScalarType >::zero())
Replace each element of the vectors in mv with alpha.
Pure virtual base class which describes the basic interface to the linear solver iteration.
static void MvTransMv(const ScalarType alpha, const MV &A, const MV &mv, Teuchos::SerialDenseMatrix< int, ScalarType > &B)
Compute a dense matrix B through the matrix-matrix multiply .
A pure virtual class for defining the status tests for the Belos iterative solvers.
Class which defines basic traits for the operator type.
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
static void MvAddMv(const ScalarType alpha, const MV &A, const ScalarType beta, const MV &B, MV &mv)
Replace mv with .
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
static Teuchos::RCP< MV > Clone(const MV &mv, const int numvecs)
Creates a new empty MV containing numvecs columns.
void iterate()
This method performs block TFQMR iterations until the status test indicates the need to stop or an er...
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.
TFQMRIterInitFailure(const std::string &what_arg)
static Teuchos::RCP< MV > CloneCopy(const MV &mv)
Creates a new MV and copies contents of mv into the new vector (deep copy).
int getNumIters() const
Get the current iteration count.
TFQMRIterInitFailure is thrown when the TFQMRIter object is unable to generate an initial iterate in ...
static ptrdiff_t GetGlobalLength(const MV &mv)
Return the number of rows in the given multivector mv.
MultiVecTraits< ScalarType, MV > MVT
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.
static void MvNorm(const MV &mv, std::vector< typename Teuchos::ScalarTraits< ScalarType >::magnitudeType > &normvec, NormType type=TwoNorm)
Compute the 2-norm of each individual vector of mv. Upon return, normvec[i] holds the value of ...
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Teuchos::RCP< const MV > D
Belos header file which uses auto-configuration information to include necessary C++ headers...
const LinearProblem< ScalarType, MV, OP > & getProblem() const
Get a constant reference to the linear problem.
Templated virtual class for providing orthogonalization/orthonormalization methods with matrix-based ...
TFQMRIterState< ScalarType, MV > getState() const
Get the current state of the linear solver.