Anasazi Version of the Day
Loading...
Searching...
No Matches
AnasaziTsqrOrthoManager.hpp
Go to the documentation of this file.
1// @HEADER
2// *****************************************************************************
3// Anasazi: Block Eigensolvers Package
4//
5// Copyright 2004 NTESS and the Anasazi contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
12
13#ifndef __AnasaziTsqrOrthoManager_hpp
14#define __AnasaziTsqrOrthoManager_hpp
15
16#include "AnasaziTsqrOrthoManagerImpl.hpp"
18
19
20namespace Anasazi {
21
43 template<class Scalar, class MV>
45 public:
46 typedef Scalar scalar_type;
47 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType magnitude_type;
50 typedef MV multivector_type;
51
52 typedef Teuchos::SerialDenseMatrix<int, Scalar> mat_type;
53 typedef Teuchos::RCP<mat_type> mat_ptr;
54
63 virtual int
64 normalizeOutOfPlace (MV& X, MV& Q, mat_ptr B) const = 0;
65
84 virtual int
86 MV& X_out,
87 Teuchos::Array<mat_ptr> C,
88 mat_ptr B,
89 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q) const = 0;
90
93 };
94
101 template<class Scalar, class MV>
103 public OrthoManager<Scalar, MV>,
104 public OutOfPlaceNormalizerMixin<Scalar, MV>,
105 public Teuchos::ParameterListAcceptor
106 {
107 public:
108 typedef Scalar scalar_type;
109 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType magnitude_type;
111 typedef MV multivector_type;
112
113 typedef Teuchos::SerialDenseMatrix<int, Scalar> mat_type;
114 typedef Teuchos::RCP<mat_type> mat_ptr;
115
116 void setParameterList (const Teuchos::RCP<Teuchos::ParameterList>& params) {
117 impl_.setParameterList (params);
118 }
119
120 Teuchos::RCP<Teuchos::ParameterList> getNonconstParameterList () {
121 return impl_.getNonconstParameterList ();
122 }
123
124 Teuchos::RCP<Teuchos::ParameterList> unsetParameterList () {
125 return impl_.unsetParameterList ();
126 }
127
135 Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const {
136 return impl_.getValidParameters();
137 }
138
148 Teuchos::RCP<const Teuchos::ParameterList> getFastParameters() const {
149 return impl_.getFastParameters();
150 }
151
168 TsqrOrthoManager (const Teuchos::RCP<Teuchos::ParameterList>& params,
169 const std::string& label = "Anasazi") :
170 impl_ (params, label)
171 {}
172
177 TsqrOrthoManager (const std::string& label) :
178 impl_ (label)
179 {}
180
182 virtual ~TsqrOrthoManager() {}
183
184 void innerProd (const MV &X, const MV& Y, mat_type& Z) const {
185 return impl_.innerProd (X, Y, Z);
186 }
187
188 void norm (const MV& X, std::vector<magnitude_type>& normVec) const {
189 return impl_.norm (X, normVec);
190 }
191
192 void
193 project (MV &X,
194 Teuchos::Array<Teuchos::RCP<const MV> > Q,
195 Teuchos::Array<Teuchos::RCP<Teuchos::SerialDenseMatrix<int,Scalar> > > C
196 = Teuchos::tuple(Teuchos::RCP< Teuchos::SerialDenseMatrix<int,Scalar> >(Teuchos::null))) const
197 {
198 return impl_.project (X, C, Q);
199 }
200
201 int
202 normalize (MV &X, mat_ptr B = Teuchos::null) const
203 {
204 return impl_.normalize (X, B);
205 }
206
207 int
209 Teuchos::Array<Teuchos::RCP<const MV> > Q,
210 Teuchos::Array<Teuchos::RCP<Teuchos::SerialDenseMatrix<int,Scalar> > > C
211 = Teuchos::tuple(Teuchos::RCP< Teuchos::SerialDenseMatrix<int,Scalar> >(Teuchos::null)),
212 Teuchos::RCP<Teuchos::SerialDenseMatrix<int,Scalar> > B = Teuchos::null) const
213 {
214 return impl_.projectAndNormalize (X, C, B, Q);
215 }
216
233 int
234 normalizeOutOfPlace (MV& X, MV& Q, mat_ptr B) const
235 {
236 return impl_.normalizeOutOfPlace (X, Q, B);
237 }
238
259 int
261 MV& X_out,
262 Teuchos::Array<mat_ptr> C,
263 mat_ptr B,
264 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q) const
265 {
266 return impl_.projectAndNormalizeOutOfPlace (X_in, X_out, C, B, Q);
267 }
268
269 magnitude_type orthonormError (const MV &X) const {
270 return impl_.orthonormError (X);
271 }
272
273 magnitude_type orthogError (const MV &X1, const MV &X2) const {
274 return impl_.orthogError (X1, X2);
275 }
276
277 private:
284 };
285
286
299 template<class Scalar, class MV, class OP>
301 public MatOrthoManager<Scalar, MV, OP>,
302 public OutOfPlaceNormalizerMixin<Scalar, MV>,
303 public Teuchos::ParameterListAcceptorDefaultBase
304 {
305 public:
306 typedef Scalar scalar_type;
307 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType magnitude_type;
311 typedef OP operator_type;
312
313 typedef Teuchos::SerialDenseMatrix<int, Scalar> mat_type;
314 typedef Teuchos::RCP<mat_type> mat_ptr;
315
316 private:
326
329 typedef TsqrOrthoManagerImpl<Scalar, MV> tsqr_type;
330
333 typedef SVQBOrthoManager<Scalar, MV, OP> svqb_type;
334
337 typedef MultiVecTraits<Scalar, MV> MVT;
338
339 public:
360 TsqrMatOrthoManager (const Teuchos::RCP<Teuchos::ParameterList>& params,
361 const std::string& label = "Belos",
362 Teuchos::RCP<const OP> Op = Teuchos::null) :
363 MatOrthoManager<Scalar, MV, OP> (Op),
364 tsqr_ (params, label),
365 pSvqb_ (Teuchos::null) // Initialized lazily
366 {}
367
376 TsqrMatOrthoManager (const std::string& label = "Belos",
377 Teuchos::RCP<const OP> Op = Teuchos::null) :
378 MatOrthoManager<Scalar, MV, OP> (Op),
379 tsqr_ (label),
380 pSvqb_ (Teuchos::null) // Initialized lazily
381 {}
382
385
393 Teuchos::RCP<const Teuchos::ParameterList> getValidParameters() const {
394 return tsqr_.getValidParameters ();
395 }
396
406 Teuchos::RCP<const Teuchos::ParameterList> getFastParameters() {
407 return tsqr_.getFastParameters ();
408 }
409
410 void setParameterList (const Teuchos::RCP<Teuchos::ParameterList>& params) {
411 tsqr_.setParameterList (params);
412 }
413
414 virtual void
415 setOp (Teuchos::RCP< const OP > Op)
416 {
417 // We override the base class' setOp() so that the
418 // SVQBOrthoManager instance gets the new op.
419 //
420 // The base_type notation helps C++ resolve the name for a
421 // member function of a templated base class.
422 base_type::setOp (Op); // base class gets a copy of the Op too
423
424 if (! Op.is_null()) {
425 ensureSvqbInit (); // Make sure the SVQB object has been initialized
426 pSvqb_->setOp (Op);
427 }
428 }
429
430 Teuchos::RCP<const OP> getOp () const {
431 // The base_type notation helps C++ resolve the name for a
432 // member function of a templated base class.
433 return base_type::getOp();
434 }
435
436 void
437 projectMat (MV &X,
438 Teuchos::Array<Teuchos::RCP<const MV> > Q,
439 Teuchos::Array<Teuchos::RCP<mat_type> > C =
440 Teuchos::tuple (Teuchos::RCP<mat_type> (Teuchos::null)),
441 Teuchos::RCP<MV> MX = Teuchos::null,
442 Teuchos::Array<Teuchos::RCP<const MV> > MQ =
443 Teuchos::tuple (Teuchos::null)) const
444 {
445 if (getOp().is_null()) {
446 // FIXME (mfh 12 Jan 2011): Do we need to check if C is null
447 // and allocate space if so?
448 tsqr_.project (X, C, Q);
449 // FIXME (mfh 20 Jul 2010) What about MX and MQ?
450 //
451 // FIXME (mfh 12 Jan 2011) Need to port MVT::Assign from Belos to Anasazi
452#if 0
453 if (! MX.is_null()) {
454 // MX gets a copy of X; M is the identity operator.
455 MVT::Assign (X, *MX);
456 }
457#endif // 0
458 } else {
459 ensureSvqbInit ();
460 pSvqb_->projectMat (X, Q, C, MX, MQ);
461 }
462 }
463
464 int
465 normalizeMat (MV &X,
466 mat_ptr B = Teuchos::null,
467 Teuchos::RCP<MV> MX = Teuchos::null) const
468 {
469 if (getOp().is_null()) {
470 // FIXME (mfh 12 Jan 2011): Do we need to check if B is null
471 // and allocate space if so?
472 const int rank = tsqr_.normalize (X, B);
473 // FIXME (mfh 20 Jul 2010) What about MX?
474 //
475 // FIXME (mfh 12 Jan 2011) Need to port MVT::Assign from Belos to Anasazi
476#if 0
477 if (! MX.is_null()) {
478 // MX gets a copy of X; M is the identity operator.
479 MVT::Assign (X, *MX);
480 }
481#endif // 0
482 return rank;
483 } else {
484 ensureSvqbInit ();
485 return pSvqb_->normalizeMat (X, B, MX);
486 }
487 }
488
489 int
491 Teuchos::Array<Teuchos::RCP<const MV> > Q,
492 Teuchos::Array<Teuchos::RCP<mat_type> > C =
493 Teuchos::tuple (Teuchos::RCP<mat_type> (Teuchos::null)),
494 Teuchos::RCP<mat_type> B = Teuchos::null,
495 Teuchos::RCP<MV> MX = Teuchos::null,
496 Teuchos::Array<Teuchos::RCP<const MV> > MQ =
497 Teuchos::tuple (Teuchos::RCP<const MV> (Teuchos::null))) const
498 {
499 if (getOp().is_null()) {
500 // FIXME (mfh 12 Jan 2011): Do we need to check if C and B
501 // are null and allocate space if so?
502 const int rank = tsqr_.projectAndNormalize (X, C, B, Q);
503 // FIXME (mfh 20 Jul 2010) What about MX and MQ?
504 // mfh 12 Jan 2011: Ignore MQ (assume Q == MQ).
505 //
506 // FIXME (mfh 12 Jan 2011) Need to port MVT::Assign from Belos to Anasazi
507#if 0
508 if (! MX.is_null()) {
509 // MX gets a copy of X; M is the identity operator.
510 MVT::Assign (X, *MX);
511 }
512#endif // 0
513 return rank;
514 } else {
515 ensureSvqbInit ();
516 return pSvqb_->projectAndNormalizeMat (X, Q, C, B, MX, MQ);
517 }
518 }
519
520 int
521 normalizeOutOfPlace (MV& X, MV& Q, mat_ptr B) const
522 {
523 if (getOp().is_null()) {
524 return tsqr_.normalizeOutOfPlace (X, Q, B);
525 } else {
526 // SVQB normalizes in place, so we have to copy.
527 ensureSvqbInit ();
528 const int rank = pSvqb_->normalize (X, B);
529 MVT::Assign (X, Q);
530 return rank;
531 }
532 }
533
534 int
536 MV& X_out,
537 Teuchos::Array<mat_ptr> C,
538 mat_ptr B,
539 Teuchos::ArrayView<Teuchos::RCP<const MV> > Q) const
540 {
541 using Teuchos::null;
542
543 if (getOp().is_null()) {
544 return tsqr_.projectAndNormalizeOutOfPlace (X_in, X_out, C, B, Q);
545 } else {
546 ensureSvqbInit ();
547 // SVQB's projectAndNormalize wants an Array, not an ArrayView.
548 // Copy into a temporary Array and copy back afterwards.
549 Teuchos::Array<Teuchos::RCP<const MV> > Q_array (Q);
550 const int rank = pSvqb_->projectAndNormalize (X_in, Q_array, C, B);
551 Q.assign (Q_array);
552 // SVQB normalizes in place, so we have to copy X_in to X_out.
553 MVT::Assign (X_in, X_out);
554 return rank;
555 }
556 }
557
558 magnitude_type
559 orthonormErrorMat (const MV &X, Teuchos::RCP<const MV> MX = Teuchos::null) const
560 {
561 if (getOp().is_null()) {
562 return tsqr_.orthonormError (X);
563 // FIXME (mfh 20 Jul 2010) What about MX?
564 } else {
565 ensureSvqbInit ();
566 return pSvqb_->orthonormErrorMat (X, MX);
567 }
568 }
569
570 magnitude_type
571 orthogErrorMat (const MV &X,
572 const MV &Y,
573 Teuchos::RCP<const MV> MX = Teuchos::null,
574 Teuchos::RCP<const MV> MY = Teuchos::null) const
575 {
576 if (getOp().is_null()) {
577 return tsqr_.orthogError (X, Y);
578 // FIXME (mfh 20 Jul 2010) What about MX and MY?
579 } else {
580 ensureSvqbInit ();
581 return pSvqb_->orthogErrorMat (X, Y, MX, MY);
582 }
583 }
584
585 private:
587 void
588 ensureSvqbInit () const
589 {
590 // NOTE (mfh 12 Oct 2011) For now, we rely on the default values
591 // of SVQB parameters.
592 if (pSvqb_.is_null()) {
593 pSvqb_ = Teuchos::rcp (new svqb_type (getOp()));
594 }
595 }
596
604 mutable tsqr_type tsqr_;
605
610 mutable Teuchos::RCP<svqb_type> pSvqb_;
611 };
612
613} // namespace Anasazi
614
615#endif // __AnasaziTsqrOrthoManager_hpp
616
Orthogonalization manager based on the SVQB technique described in "A Block Orthogonalization Procedu...
Anasazi's templated virtual class for providing routines for orthogonalization and orthonormalization...
virtual Teuchos::RCP< const OP > getOp() const
Get operator used for inner product.
virtual void setOp(Teuchos::RCP< const OP > Op)
Set operator used for inner product.
Traits class which defines basic operations on multivectors.
Anasazi's templated virtual class for providing routines for orthogonalization and orthonormalization...
Mixin for out-of-place orthogonalization.
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.
virtual ~OutOfPlaceNormalizerMixin()
Trivial virtual destructor, to silence compiler warnings.
virtual int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const =0
Normalize X into Q*B.
MV multivector_type
Multivector type with which this class was specialized.
An implementation of the Anasazi::MatOrthoManager that performs orthogonalization using the SVQB iter...
MatOrthoManager subclass using TSQR or SVQB.
int projectAndNormalizeMat(MV &X, Teuchos::Array< Teuchos::RCP< const MV > > Q, Teuchos::Array< Teuchos::RCP< mat_type > > C=Teuchos::tuple(Teuchos::RCP< mat_type >(Teuchos::null)), Teuchos::RCP< mat_type > B=Teuchos::null, Teuchos::RCP< MV > MX=Teuchos::null, Teuchos::Array< Teuchos::RCP< const MV > > MQ=Teuchos::tuple(Teuchos::RCP< const MV >(Teuchos::null))) const
Provides matrix-based projection/orthonormalization method.
Teuchos::RCP< const Teuchos::ParameterList > getFastParameters()
Get "fast" parameters for TsqrMatOrthoManager.
OP operator_type
Operator type with which this class was specialized.
virtual void setOp(Teuchos::RCP< const OP > Op)
Set operator used for inner product.
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.
int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const
Normalize X into Q*B.
TsqrMatOrthoManager(const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null)
Constructor (that sets default parameters).
TsqrMatOrthoManager(const Teuchos::RCP< Teuchos::ParameterList > &params, const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null)
Constructor (that sets user-specified parameters).
magnitude_type orthonormErrorMat(const MV &X, Teuchos::RCP< const MV > MX=Teuchos::null) const
This method computes the error in orthonormality of a multivector.
magnitude_type orthogErrorMat(const MV &X, const MV &Y, Teuchos::RCP< const MV > MX=Teuchos::null, Teuchos::RCP< const MV > MY=Teuchos::null) const
This method computes the error in orthogonality of two multivectors.
Teuchos::RCP< const OP > getOp() const
Get operator used for inner product.
MV multivector_type
Multivector type with which this class was specialized.
void projectMat(MV &X, Teuchos::Array< Teuchos::RCP< const MV > > Q, Teuchos::Array< Teuchos::RCP< mat_type > > C=Teuchos::tuple(Teuchos::RCP< mat_type >(Teuchos::null)), Teuchos::RCP< MV > MX=Teuchos::null, Teuchos::Array< Teuchos::RCP< const MV > > MQ=Teuchos::tuple(Teuchos::null)) const
Provides matrix-based projection method.
virtual ~TsqrMatOrthoManager()
Destructor (declared virtual for memory safety of derived classes).
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Get default parameters for TsqrMatOrthoManager.
TSQR-based OrthoManager subclass implementation.
TSQR-based OrthoManager subclass.
magnitude_type orthogError(const MV &X1, const MV &X2) const
This method computes the error in orthogonality of two multivectors.
virtual ~TsqrOrthoManager()
Destructor, declared virtual for safe inheritance.
int projectAndNormalize(MV &X, Teuchos::Array< Teuchos::RCP< const MV > > Q, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, Scalar > > > C=Teuchos::tuple(Teuchos::RCP< Teuchos::SerialDenseMatrix< int, Scalar > >(Teuchos::null)), Teuchos::RCP< Teuchos::SerialDenseMatrix< int, Scalar > > B=Teuchos::null) const
Given a set of bases Q[i] and a multivector X, this method computes an orthonormal basis for .
Teuchos::RCP< const Teuchos::ParameterList > getFastParameters() const
Get "fast" parameters for TsqrOrthoManager.
magnitude_type orthonormError(const MV &X) const
This method computes the error in orthonormality of a multivector.
TsqrOrthoManager(const Teuchos::RCP< Teuchos::ParameterList > &params, const std::string &label="Anasazi")
Constructor (that sets user-specified parameters).
TsqrOrthoManager(const std::string &label)
Constructor (that sets default parameters).
int normalizeOutOfPlace(MV &X, MV &Q, mat_ptr B) const
Normalize X into Q*B, overwriting X with invalid values.
void innerProd(const MV &X, const MV &Y, mat_type &Z) const
Provides the inner product defining the orthogonality concepts.
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< const Teuchos::ParameterList > getValidParameters() const
Default valid parameter list.
void project(MV &X, Teuchos::Array< Teuchos::RCP< const MV > > Q, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, Scalar > > > C=Teuchos::tuple(Teuchos::RCP< Teuchos::SerialDenseMatrix< int, Scalar > >(Teuchos::null))) const
Given a list of mutually orthogonal and internally orthonormal bases Q, this method projects a multiv...
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package.