45#ifndef __BelosTsqrOrthoManager_hpp
46#define __BelosTsqrOrthoManager_hpp
78template<
class Scalar,
class MV>
82 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType
magnitude_type;
87 typedef Teuchos::SerialDenseMatrix<int, Scalar>
mat_type;
122 Teuchos::Array<mat_ptr> C,
124 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const = 0;
136template<
class Scalar,
class MV>
147 typedef Teuchos::SerialDenseMatrix<int, Scalar>
mat_type;
151 impl_.setParameterList (params);
155 return impl_.getNonconstParameterList ();
159 return impl_.unsetParameterList ();
170 return impl_.getValidParameters();
183 return impl_.getFastParameters();
202 const std::string& label =
"Belos") :
203 impl_ (params, label)
241 impl_.setReorthogonalizationCallback (callback);
245 return impl_.innerProd (X, Y, Z);
248 void norm (
const MV& X, std::vector<magnitude_type>& normVec)
const {
249 return impl_.norm (X, normVec);
254 Teuchos::Array<mat_ptr> C,
255 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
257 return impl_.project (X, C, Q);
263 return impl_.normalize (X, B);
269 Teuchos::Array<mat_ptr> C,
271 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
273 return impl_.projectAndNormalize (X, C, B, Q);
296 return impl_.normalizeOutOfPlace (X, Q, B);
322 Teuchos::Array<mat_ptr> C,
324 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
326 return impl_.projectAndNormalizeOutOfPlace (X_in, X_out, C, B, Q);
330 return impl_.orthonormError (X);
334 return impl_.orthogError (X1, X2);
345 impl_.setLabel (label);
376template<
class Scalar,
class MV,
class OP>
389 typedef Teuchos::SerialDenseMatrix<int, Scalar>
mat_type;
437 const std::string& label =
"Belos",
438 Teuchos::RCP<const OP> Op = Teuchos::null) :
440 tsqr_ (params, label),
453 Teuchos::RCP<const OP> Op = Teuchos::null) :
502 if (! Op.is_null()) {
508 Teuchos::RCP<const OP>
getOp ()
const {
517 Teuchos::Array<mat_ptr> C,
518 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
520 if (
getOp().is_null()) {
522 if (! MX.is_null()) {
528 pDgks_->project (X, MX, C, Q);
534 Teuchos::Array<mat_ptr> C,
535 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
537 project (X, Teuchos::null, C, Q);
543 if (
getOp().is_null()) {
545 if (! MX.is_null()) {
552 return pDgks_->normalize (X, MX, B);
570 Teuchos::Array<mat_ptr> C,
572 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
574 if (
getOp().is_null()) {
576 if (! MX.is_null()) {
583 return pDgks_->projectAndNormalize (X, MX, C, B, Q);
591 if (
getOp().is_null()) {
596 const int rank =
pDgks_->normalize (X, B);
605 Teuchos::Array<mat_ptr> C,
607 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q)
const
611 if (
getOp().is_null()) {
616 const int rank =
pDgks_->projectAndNormalize (X_in, null, C, B, Q);
625 if (
getOp().is_null()) {
629 return pDgks_->orthonormError (X, MX);
643 Teuchos::RCP<const MV> MX1,
646 if (
getOp ().is_null ()) {
651 return pDgks_->orthogError (X1, MX1, X2);
662 if (!
pDgks_.is_null ()) {
Classical Gram-Schmidt (with DGKS correction) implementation of the Belos::OrthoManager class.
Orthogonalization manager back end based on Tall Skinny QR (TSQR)
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
Belos's templated virtual class for providing routines for orthogonalization and orthonormzalition of...
void setOp(Teuchos::RCP< const OP > Op)
Set operator.
Teuchos::RCP< const OP > getOp() const
Get operator.
Traits class which defines basic operations on multivectors.
static void Assign(const MV &A, MV &mv)
mv := A
Belos's templated virtual class for providing routines for orthogonalization and orthonormzalition of...
Mixin for out-of-place orthogonalization.
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitude_type
MV multivector_type
Multivector type with which this class was specialized.
virtual int projectAndNormalizeOutOfPlace(MV &X_in, MV &X_out, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const =0
Project and normalize X_in into X_out.
Teuchos::RCP< mat_type > mat_ptr
virtual int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const =0
Normalize X into Q*B.
Teuchos::SerialDenseMatrix< int, Scalar > mat_type
virtual ~OutOfPlaceNormalizerMixin()
Trivial virtual destructor, to silence compiler warnings.
Interface of callback invoked by TsqrOrthoManager on reorthogonalization.
MatOrthoManager subclass using TSQR or DGKS.
void ensureDgksInit() const
Ensure that the DGKSOrthoManager instance is initialized.
int normalize(MV &X, Teuchos::RCP< MV > MX, mat_ptr B) const
virtual ~TsqrMatOrthoManager()
Destructor (declared virtual for memory safety of derived classes).
MultiVecTraits< Scalar, MV > MVT
Traits class for the multivector type.
TsqrOrthoManagerImpl< Scalar, MV > tsqr_type
Implementation of TSQR-based orthogonalization.
magnitude_type orthogError(const MV &X1, Teuchos::RCP< const MV > MX1, const MV &X2) const
This method computes the error in orthogonality of two multivectors. The method has the option of exp...
Teuchos::RCP< const Teuchos::ParameterList > getFastParameters()
Get "fast" parameters for TsqrMatOrthoManager.
int normalize(MV &X, mat_ptr B) const
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
void project(MV &X, Teuchos::Array< mat_ptr > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Teuchos::SerialDenseMatrix< int, Scalar > mat_type
TsqrMatOrthoManager(const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null)
Constructor (that sets default parameters).
magnitude_type orthonormError(const MV &X) const
This method computes the error in orthonormality of a multivector.
int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const
Normalize X into Q*B.
tsqr_type tsqr_
TSQR + BGS orthogonalization manager implementation.
void setLabel(const std::string &label)
This method sets the label used by the timers in the orthogonalization manager.
Teuchos::RCP< const OP > getOp() const
MV multivector_type
Multivector type with which this class was specialized.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Get default parameters for TsqrMatOrthoManager.
MatOrthoManager< Scalar, MV, OP > base_type
This will be used to help C++ resolve getOp().
magnitude_type orthogError(const MV &X1, const MV &X2) const
This method computes the error in orthogonality of two multivectors. This method.
DGKSOrthoManager< Scalar, MV, OP > dgks_type
Implementation of DGKS-based orthogonalization.
void setOp(Teuchos::RCP< const OP > Op)
const std::string & getLabel() const
This method returns the label being used by the timers in the orthogonalization manager.
magnitude_type orthonormError(const MV &X, Teuchos::RCP< const MV > MX) const
This method computes the error in orthonormality of a multivector. The method has the option of explo...
Teuchos::RCP< dgks_type > pDgks_
DGKS orthogonalization manager.
TsqrMatOrthoManager(const Teuchos::RCP< Teuchos::ParameterList > ¶ms, const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null)
Constructor (that sets user-specified parameters).
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitude_type
void project(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< mat_ptr > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
int projectAndNormalizeOutOfPlace(MV &X_in, MV &X_out, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Project and normalize X_in into X_out.
Teuchos::RCP< mat_type > mat_ptr
OP operator_type
Operator type with which this class was specialized.
virtual int projectAndNormalizeWithMxImpl(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
TSQR-based OrthoManager subclass implementation.
int projectAndNormalizeOutOfPlace(MV &X_in, MV &X_out, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q)
Project and normalize X_in into X_out; overwrite X_in.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Set parameters from the given parameter list.
int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B)
Normalize X into Q*B, overwriting X.
void setLabel(const std::string &label)
Set the label for timers.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Default valid parameter list.
int projectAndNormalize(MV &X, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q)
Project X against Q and normalize X.
magnitude_type orthogError(const MV &X1, const MV &X2) const
Return the Frobenius norm of the inner product of X1 with itself.
void project(MV &X, Teuchos::Array< mat_ptr > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q)
Compute and .
int normalize(MV &X, mat_ptr B)
Orthogonalize the columns of X in place.
Teuchos::RCP< const Teuchos::ParameterList > getFastParameters()
Get "fast" parameters for TsqrOrthoManagerImpl.
const std::string & getLabel() const
Get the label for timers (if timers are enabled).
magnitude_type orthonormError(const MV &X) const
Return .
TSQR-based OrthoManager subclass.
TsqrOrthoManager(const std::string &label)
Constructor (that sets default parameters).
magnitude_type orthonormError(const MV &X) const
This method computes the error in orthonormality of a multivector.
void project(MV &X, Teuchos::Array< mat_ptr > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
void norm(const MV &X, std::vector< magnitude_type > &normVec) const
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitude_type
const std::string & getLabel() const
This method returns the label being used by the timers in the orthogonalization manager.
Teuchos::RCP< const Teuchos::ParameterList > getFastParameters() const
Get "fast" parameters for TsqrOrthoManager.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
void innerProd(const MV &X, const MV &Y, mat_type &Z) const
Provides the inner product defining the orthogonality concepts.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Default valid parameter list.
virtual int projectAndNormalizeImpl(MV &X, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
magnitude_type orthogError(const MV &X1, const MV &X2) const
This method computes the error in orthogonality of two multivectors.
int normalize(MV &X, mat_ptr B) const
int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const
Normalize X into Q*B, overwriting X with invalid values.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
TsqrOrthoManager(const Teuchos::RCP< Teuchos::ParameterList > ¶ms, const std::string &label="Belos")
Constructor (that sets user-specified parameters).
TsqrOrthoManagerImpl< Scalar, MV > impl_
The implementation of TSQR.
std::string label_
Label for timers (if timers are enabled at compile time).
Teuchos::SerialDenseMatrix< int, Scalar > mat_type
virtual ~TsqrOrthoManager()
Destructor, declared virtual for safe inheritance.
int projectAndNormalizeOutOfPlace(MV &X_in, MV &X_out, Teuchos::Array< mat_ptr > C, mat_ptr B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Project and normalize X_in into X_out; overwrite X_in.
Teuchos::RCP< mat_type > mat_ptr
void setReorthogonalizationCallback(const Teuchos::RCP< ReorthogonalizationCallback< Scalar > > &callback)
Set callback to be invoked on reorthogonalization.
void setLabel(const std::string &label)
Set the label for (the timers for) this orthogonalization manager, and create new timers if the label...