47 #ifndef BELOS_ICGS_ORTHOMANAGER_HPP
48 #define BELOS_ICGS_ORTHOMANAGER_HPP
65 #ifdef BELOS_TEUCHOS_TIME_MONITOR
67 #endif // BELOS_TEUCHOS_TIME_MONITOR
72 template<
class ScalarType,
class MV,
class OP>
76 template<
class ScalarType,
class MV,
class OP>
79 template<
class ScalarType,
class MV,
class OP>
106 #ifdef BELOS_TEUCHOS_TIME_MONITOR
107 std::stringstream ss;
110 std::string orthoLabel = ss.str() +
": Orthogonalization";
113 std::string updateLabel = ss.str() +
": Ortho (Update)";
116 std::string normLabel = ss.str() +
": Ortho (Norm)";
119 std::string ipLabel = ss.str() +
": Ortho (Inner Product)";
126 const std::string& label =
"Belos",
136 #ifdef BELOS_TEUCHOS_TIME_MONITOR
137 std::stringstream ss;
140 std::string orthoLabel = ss.str() +
": Orthogonalization";
143 std::string updateLabel = ss.str() +
": Ortho (Update)";
146 std::string normLabel = ss.str() +
": Ortho (Norm)";
149 std::string ipLabel = ss.str() +
": Ortho (Inner Product)";
166 using Teuchos::parameterList;
170 RCP<ParameterList> params;
172 params = parameterList (*defaultParams);
187 int maxNumOrthogPasses;
192 maxNumOrthogPasses = params->
get<
int> (
"maxNumOrthogPasses");
193 }
catch (InvalidParameterName&) {
194 maxNumOrthogPasses = defaultParams->get<
int> (
"maxNumOrthogPasses");
195 params->set (
"maxNumOrthogPasses", maxNumOrthogPasses);
207 }
catch (InvalidParameterName&) {
212 params->remove (
"depTol");
213 }
catch (InvalidParameterName&) {
216 params->set (
"blkTol", blkTol);
221 }
catch (InvalidParameterName&) {
223 params->set (
"singTol", singTol);
237 defaultParams_ = Belos::getICGSDefaultParameters<ScalarType, MV, OP>();
254 using Teuchos::parameterList;
259 RCP<ParameterList> params = parameterList (*defaultParams);
280 params->
set (
"blkTol", blk_tol);
294 params->
set (
"singTol", sing_tol);
479 void setLabel(
const std::string& label);
517 #ifdef BELOS_TEUCHOS_TIME_MONITOR
519 #endif // BELOS_TEUCHOS_TIME_MONITOR
527 bool completeBasis,
int howMany = -1 )
const;
559 template<
class ScalarType,
class MV,
class OP>
562 template<
class ScalarType,
class MV,
class OP>
568 template<
class ScalarType,
class MV,
class OP>
573 template<
class ScalarType,
class MV,
class OP>
576 template<
class ScalarType,
class MV,
class OP>
581 template<
class ScalarType,
class MV,
class OP>
588 template<
class ScalarType,
class MV,
class OP>
591 if (label != label_) {
593 #ifdef BELOS_TEUCHOS_TIME_MONITOR
594 std::stringstream ss;
595 ss << label_ +
": ICGS[" << max_ortho_steps_ <<
"]";
597 std::string orthoLabel = ss.str() +
": Orthogonalization";
600 std::string updateLabel = ss.str() +
": Ortho (Update)";
603 std::string normLabel = ss.str() +
": Ortho (Norm)";
606 std::string ipLabel = ss.str() +
": Ortho (Inner Product)";
614 template<
class ScalarType,
class MV,
class OP>
617 const ScalarType ONE = SCT::one();
618 int rank = MVT::GetNumberVecs(X);
621 for (
int i=0; i<rank; i++) {
629 template<
class ScalarType,
class MV,
class OP>
632 int r1 = MVT::GetNumberVecs(X1);
633 int r2 = MVT::GetNumberVecs(X2);
641 template<
class ScalarType,
class MV,
class OP>
657 typedef typename Array< RCP< const MV > >::size_type size_type;
659 #ifdef BELOS_TEUCHOS_TIME_MONITOR
663 ScalarType ONE = SCT::one();
667 int xc = MVT::GetNumberVecs( X );
668 ptrdiff_t xr = MVT::GetGlobalLength( X );
675 B =
rcp (
new serial_dense_matrix_type (xc, xc));
685 for (size_type k = 0; k < nq; ++k)
687 const int numRows = MVT::GetNumberVecs (*Q[k]);
688 const int numCols = xc;
691 C[k] =
rcp (
new serial_dense_matrix_type (numRows, numCols));
692 else if (
C[k]->numRows() != numRows ||
C[k]->numCols() != numCols)
694 int err =
C[k]->reshape (numRows, numCols);
696 "IMGS orthogonalization: failed to reshape "
697 "C[" << k <<
"] (the array of block "
698 "coefficients resulting from projecting X "
699 "against Q[1:" << nq <<
"]).");
707 MX = MVT::Clone(X,MVT::GetNumberVecs(X));
708 OPT::Apply(*(this->_Op),X,*MX);
716 int mxc = MVT::GetNumberVecs( *MX );
717 ptrdiff_t mxr = MVT::GetGlobalLength( *MX );
720 TEUCHOS_TEST_FOR_EXCEPTION( xc == 0 || xr == 0, std::invalid_argument,
"Belos::ICGSOrthoManager::projectAndNormalize(): X must be non-empty" );
723 for (
int i=0; i<nq; i++) {
724 numbas += MVT::GetNumberVecs( *Q[i] );
729 "Belos::ICGSOrthoManager::projectAndNormalize(): Size of X must be consistant with size of B" );
732 "Belos::ICGSOrthoManager::projectAndNormalize(): MVT returned negative dimensions for X,MX" );
734 TEUCHOS_TEST_FOR_EXCEPTION( xc!=mxc || xr!=mxr, std::invalid_argument,
735 "Belos::ICGSOrthoManager::projectAndNormalize(): Size of X must be consistant with size of MX" );
741 bool dep_flg =
false;
747 dep_flg = blkOrtho1( X, MX,
C, Q );
753 std::vector<ScalarType> diag(xc);
755 #ifdef BELOS_TEUCHOS_TIME_MONITOR
758 MVT::MvDot( X, *MX, diag );
760 (*B)(0,0) = SCT::squareroot(SCT::magnitude(diag[0]));
762 if (SCT::magnitude((*
B)(0,0)) > ZERO) {
764 MVT::MvScale( X, ONE/(*
B)(0,0) );
767 MVT::MvScale( *MX, ONE/(*
B)(0,0) );
775 tmpX = MVT::CloneCopy(X);
777 tmpMX = MVT::CloneCopy(*MX);
781 dep_flg = blkOrtho( X, MX,
C, Q );
787 rank = blkOrthoSing( *tmpX, tmpMX,
C,
B, Q );
790 MVT::Assign( *tmpX, X );
792 MVT::Assign( *tmpMX, *MX );
797 rank = findBasis( X, MX,
B,
false );
802 rank = blkOrthoSing( *tmpX, tmpMX,
C,
B, Q );
805 MVT::Assign( *tmpX, X );
807 MVT::Assign( *tmpMX, *MX );
814 TEUCHOS_TEST_FOR_EXCEPTION( rank > xc || rank < 0, std::logic_error,
815 "Belos::ICGSOrthoManager::projectAndNormalize(): Debug error in rank variable." );
825 template<
class ScalarType,
class MV,
class OP>
830 #ifdef BELOS_TEUCHOS_TIME_MONITOR
835 return findBasis(X, MX,
B,
true);
842 template<
class ScalarType,
class MV,
class OP>
862 #ifdef BELOS_TEUCHOS_TIME_MONITOR
866 int xc = MVT::GetNumberVecs( X );
867 ptrdiff_t xr = MVT::GetGlobalLength( X );
869 std::vector<int> qcs(nq);
871 if (nq == 0 || xc == 0 || xr == 0) {
874 ptrdiff_t qr = MVT::GetGlobalLength ( *Q[0] );
885 MX = MVT::Clone(X,MVT::GetNumberVecs(X));
886 OPT::Apply(*(this->_Op),X,*MX);
893 int mxc = MVT::GetNumberVecs( *MX );
894 ptrdiff_t mxr = MVT::GetGlobalLength( *MX );
898 "Belos::ICGSOrthoManager::project(): MVT returned negative dimensions for X,MX" );
901 "Belos::ICGSOrthoManager::project(): Size of X not consistant with MX,Q" );
904 for (
int i=0; i<nq; i++) {
906 "Belos::ICGSOrthoManager::project(): Q lengths not mutually consistant" );
907 qcs[i] = MVT::GetNumberVecs( *Q[i] );
909 "Belos::ICGSOrthoManager::project(): Q has less rows than columns" );
917 "Belos::ICGSOrthoManager::project(): Size of Q not consistant with size of C" );
922 blkOrtho( X, MX,
C, Q );
929 template<
class ScalarType,
class MV,
class OP>
952 const ScalarType ONE = SCT::one ();
955 const int xc = MVT::GetNumberVecs (X);
956 const ptrdiff_t xr = MVT::GetGlobalLength (X);
971 MX = MVT::Clone(X,xc);
972 OPT::Apply(*(this->_Op),X,*MX);
983 const int mxc = (this->_hasOp) ? MVT::GetNumberVecs( *MX ) : xc;
984 const ptrdiff_t mxr = (this->_hasOp) ? MVT::GetGlobalLength( *MX ) : xr;
988 "Belos::ICGSOrthoManager::findBasis(): X must be non-empty" );
990 "Belos::ICGSOrthoManager::findBasis(): Size of X not consistant with size of B" );
992 "Belos::ICGSOrthoManager::findBasis(): Size of X not consistant with size of MX" );
993 TEUCHOS_TEST_FOR_EXCEPTION( as<ptrdiff_t> (xc) > xr, std::invalid_argument,
994 "Belos::ICGSOrthoManager::findBasis(): Size of X not feasible for normalization" );
995 TEUCHOS_TEST_FOR_EXCEPTION( howMany < 0 || howMany > xc, std::invalid_argument,
996 "Belos::ICGSOrthoManager::findBasis(): Invalid howMany parameter" );
1001 int xstart = xc - howMany;
1003 for (
int j = xstart; j < xc; j++) {
1009 bool addVec =
false;
1012 std::vector<int> index(1);
1018 MXj = MVT::CloneViewNonConst( *MX, index );
1026 std::vector<int> prev_idx( numX );
1031 for (
int i=0; i<numX; i++) {
1034 prevX = MVT::CloneView( X, prev_idx );
1036 prevMX = MVT::CloneView( *MX, prev_idx );
1039 oldMXj = MVT::CloneCopy( *MXj );
1044 std::vector<ScalarType> oldDot( 1 ), newDot( 1 );
1049 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1052 MVT::MvDot( *Xj, *MXj, oldDot );
1055 TEUCHOS_TEST_FOR_EXCEPTION( SCT::real(oldDot[0]) < ZERO,
OrthoError,
1056 "Belos::ICGSOrthoManager::findBasis(): Negative definiteness discovered in inner product" );
1062 for (
int i=0; i<max_ortho_steps_; ++i) {
1066 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1075 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1078 MVT::MvTimesMatAddMv( -ONE, *prevX, P2, ONE, *Xj );
1086 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1089 MVT::MvTimesMatAddMv( -ONE, *prevMX, P2, ONE, *MXj );
1103 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1106 MVT::MvDot( *Xj, *oldMXj, newDot );
1109 newDot[0] = oldDot[0];
1113 if (completeBasis) {
1117 if ( SCT::magnitude(newDot[0]) < SCT::magnitude(sing_tol_*oldDot[0]) ) {
1122 std::cout <<
"Belos::ICGSOrthoManager::findBasis() --> Random for column " << numX << std::endl;
1127 MVT::MvRandom( *tempXj );
1129 tempMXj = MVT::Clone( X, 1 );
1130 OPT::Apply( *(this->_Op), *tempXj, *tempMXj );
1136 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1139 MVT::MvDot( *tempXj, *tempMXj, oldDot );
1142 for (
int num_orth=0; num_orth<max_ortho_steps_; num_orth++){
1144 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1150 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1153 MVT::MvTimesMatAddMv( -ONE, *prevX, product, ONE, *tempXj );
1156 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1159 MVT::MvTimesMatAddMv( -ONE, *prevMX, product, ONE, *tempMXj );
1164 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1167 MVT::MvDot( *tempXj, *tempMXj, newDot );
1170 if ( SCT::magnitude(newDot[0]) >= SCT::magnitude(oldDot[0]*sing_tol_) ) {
1172 MVT::Assign( *tempXj, *Xj );
1174 MVT::Assign( *tempMXj, *MXj );
1186 if ( SCT::magnitude(newDot[0]) < SCT::magnitude(oldDot[0]*blk_tol_) ) {
1194 ScalarType diag = SCT::squareroot(SCT::magnitude(newDot[0]));
1195 if (SCT::magnitude(diag) > ZERO) {
1196 MVT::MvScale( *Xj, ONE/diag );
1199 MVT::MvScale( *MXj, ONE/diag );
1213 for (
int i=0; i<numX; i++) {
1214 (*B)(i,j) = product(i,0);
1225 template<
class ScalarType,
class MV,
class OP>
1232 int xc = MVT::GetNumberVecs( X );
1233 const ScalarType ONE = SCT::one();
1235 std::vector<int> qcs( nq );
1236 for (
int i=0; i<nq; i++) {
1237 qcs[i] = MVT::GetNumberVecs( *Q[i] );
1244 for (
int i=0; i<nq; i++) {
1247 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1254 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1257 MVT::MvTimesMatAddMv( -ONE, *Q[i], *
C[i], ONE, X );
1263 OPT::Apply( *(this->_Op), X, *MX);
1267 MQ[i] = MVT::Clone( *Q[i], qcs[i] );
1268 OPT::Apply( *(this->_Op), *Q[i], *MQ[i] );
1270 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1273 MVT::MvTimesMatAddMv( -ONE, *MQ[i], *
C[i], ONE, *MX );
1280 for (
int j = 1; j < max_ortho_steps_; ++j) {
1282 for (
int i=0; i<nq; i++) {
1287 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1294 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1297 MVT::MvTimesMatAddMv( -ONE, *Q[i], C2, ONE, X );
1303 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1307 MVT::MvTimesMatAddMv( -ONE, *MQ[i], C2, ONE, *MX );
1309 else if (xc <= qcs[i]) {
1311 OPT::Apply( *(this->_Op), X, *MX);
1322 template<
class ScalarType,
class MV,
class OP>
1329 int xc = MVT::GetNumberVecs( X );
1330 bool dep_flg =
false;
1331 const ScalarType ONE = SCT::one();
1333 std::vector<int> qcs( nq );
1334 for (
int i=0; i<nq; i++) {
1335 qcs[i] = MVT::GetNumberVecs( *Q[i] );
1341 std::vector<ScalarType> oldDot( xc );
1343 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1346 MVT::MvDot( X, *MX, oldDot );
1351 for (
int i=0; i<nq; i++) {
1354 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1361 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1364 MVT::MvTimesMatAddMv( -ONE, *Q[i], *
C[i], ONE, X );
1369 OPT::Apply( *(this->_Op), X, *MX);
1373 MQ[i] = MVT::Clone( *Q[i], qcs[i] );
1374 OPT::Apply( *(this->_Op), *Q[i], *MQ[i] );
1376 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1379 MVT::MvTimesMatAddMv( -ONE, *MQ[i], *
C[i], ONE, *MX );
1386 for (
int j = 1; j < max_ortho_steps_; ++j) {
1388 for (
int i=0; i<nq; i++) {
1393 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1400 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1403 MVT::MvTimesMatAddMv( -ONE, *Q[i], C2, ONE, X );
1409 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1413 MVT::MvTimesMatAddMv( -ONE, *MQ[i], C2, ONE, *MX );
1415 else if (xc <= qcs[i]) {
1417 OPT::Apply( *(this->_Op), X, *MX);
1424 std::vector<ScalarType> newDot(xc);
1426 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1429 MVT::MvDot( X, *MX, newDot );
1433 for (
int i=0; i<xc; i++){
1434 if (SCT::magnitude(newDot[i]) < SCT::magnitude(oldDot[i] * blk_tol_)) {
1443 template<
class ScalarType,
class MV,
class OP>
1452 const ScalarType ONE = SCT::one();
1453 const ScalarType ZERO = SCT::zero();
1456 int xc = MVT::GetNumberVecs( X );
1457 std::vector<int> indX( 1 );
1458 std::vector<ScalarType> oldDot( 1 ), newDot( 1 );
1460 std::vector<int> qcs( nq );
1461 for (
int i=0; i<nq; i++) {
1462 qcs[i] = MVT::GetNumberVecs( *Q[i] );
1471 for (
int j=0; j<xc; j++) {
1473 bool dep_flg =
false;
1477 std::vector<int> index( j );
1478 for (
int ind=0; ind<j; ind++) {
1481 lastQ = MVT::CloneView( X, index );
1486 qcs.push_back( MVT::GetNumberVecs( *lastQ ) );
1491 Xj = MVT::CloneViewNonConst( X, indX );
1493 MXj = MVT::CloneViewNonConst( *MX, indX );
1501 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1504 MVT::MvDot( *Xj, *MXj, oldDot );
1509 for (
int i=0; i<Q.
size(); i++) {
1516 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1522 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1526 MVT::MvTimesMatAddMv( -ONE, *Q[i], tempC, ONE, *Xj );
1531 OPT::Apply( *(this->_Op), *Xj, *MXj);
1535 MQ[i] = MVT::Clone( *Q[i], qcs[i] );
1536 OPT::Apply( *(this->_Op), *Q[i], *MQ[i] );
1538 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1541 MVT::MvTimesMatAddMv( -ONE, *MQ[i], tempC, ONE, *MXj );
1548 for (
int num_ortho_steps=1; num_ortho_steps < max_ortho_steps_; ++num_ortho_steps) {
1550 for (
int i=0; i<Q.
size(); i++) {
1556 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1563 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1566 MVT::MvTimesMatAddMv( -ONE, *Q[i], C2, ONE, *Xj );
1573 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1576 MVT::MvTimesMatAddMv( -ONE, *MQ[i], C2, ONE, *MXj );
1578 else if (xc <= qcs[i]) {
1580 OPT::Apply( *(this->_Op), *Xj, *MXj);
1589 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1592 MVT::MvDot( *Xj, *MXj, newDot );
1596 if (SCT::magnitude(newDot[0]) < SCT::magnitude(oldDot[0]*sing_tol_)) {
1602 ScalarType diag = SCT::squareroot(SCT::magnitude(newDot[0]));
1604 MVT::MvScale( *Xj, ONE/diag );
1607 MVT::MvScale( *MXj, ONE/diag );
1617 MVT::MvRandom( *tempXj );
1619 tempMXj = MVT::Clone( X, 1 );
1620 OPT::Apply( *(this->_Op), *tempXj, *tempMXj );
1626 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1629 MVT::MvDot( *tempXj, *tempMXj, oldDot );
1632 for (
int num_orth=0; num_orth<max_ortho_steps_; num_orth++) {
1634 for (
int i=0; i<Q.
size(); i++) {
1639 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1645 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1648 MVT::MvTimesMatAddMv( -ONE, *Q[i], product, ONE, *tempXj );
1654 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1658 MVT::MvTimesMatAddMv( -ONE, *MQ[i], product, ONE, *tempMXj );
1660 else if (xc <= qcs[i]) {
1662 OPT::Apply( *(this->_Op), *tempXj, *tempMXj);
1671 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1674 MVT::MvDot( *tempXj, *tempMXj, newDot );
1678 if ( SCT::magnitude(newDot[0]) >= SCT::magnitude(oldDot[0]*sing_tol_) ) {
1679 ScalarType diag = SCT::squareroot(SCT::magnitude(newDot[0]));
1685 MVT::MvAddMv( ONE/diag, *tempXj, ZERO, *tempXj, *Xj );
1687 MVT::MvAddMv( ONE/diag, *tempMXj, ZERO, *tempMXj, *MXj );
1707 template<
class ScalarType,
class MV,
class OP>
1711 using Teuchos::parameterList;
1714 RCP<ParameterList> params = parameterList (
"ICGS");
1719 "Maximum number of orthogonalization passes (includes the "
1720 "first). Default is 2, since \"twice is enough\" for Krylov "
1723 "Block reorthogonalization threshold.");
1725 "Singular block detection threshold.");
1730 template<
class ScalarType,
class MV,
class OP>
1736 RCP<ParameterList> params = getICGSDefaultParameters<ScalarType, MV, OP>();
1738 params->set (
"maxNumOrthogPasses",
1740 params->set (
"blkTol",
1742 params->set (
"singTol",
1750 #endif // BELOS_ICGS_ORTHOMANAGER_HPP
void setBlkTol(const MagnitudeType blk_tol)
Set parameter for block re-orthogonalization threshhold.
const std::string & getLabel() const
This method returns the label being used by the timers in the orthogonalization manager.
void setParameterList(const Teuchos::RCP< Teuchos::ParameterList > &plist)
bool is_null(const boost::shared_ptr< T > &p)
static const MagnitudeType sing_tol_default_
Singular block detection threshold (default).
ICGSOrthoManager(const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null, const int max_ortho_steps=max_ortho_steps_default_, const MagnitudeType blk_tol=blk_tol_default_, const MagnitudeType sing_tol=sing_tol_default_)
Constructor specifying re-orthogonalization tolerance.
static const int max_ortho_steps_default_
Max number of (re)orthogonalization steps, including the first (default).
T & get(ParameterList &l, const std::string &name)
Teuchos::RCP< Teuchos::ParameterList > getICGSDefaultParameters()
"Default" parameters for robustness and accuracy.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
bool is_null(const std::shared_ptr< T > &p)
ICGSOrthoManager(const Teuchos::RCP< Teuchos::ParameterList > &plist, const std::string &label="Belos", Teuchos::RCP< const OP > Op=Teuchos::null)
Constructor that takes a list of parameters.
MultiVecTraits< ScalarType, MV > MVT
MagnitudeType getBlkTol() const
Return parameter for block re-orthogonalization threshhold.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Declaration of basic traits for the multivector type.
void setLabel(const std::string &label)
This method sets the label used by the timers in the orthogonalization manager.
int blkOrthoSing(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B, Teuchos::ArrayView< Teuchos::RCP< const MV > > QQ) const
Project X against QQ and normalize X, one vector at a time.
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
Teuchos::ScalarTraits< ScalarType >::magnitudeType orthogError(const MV &X1, const MV &X2) const
This method computes the error in orthogonality of two multivectors, measured as the Frobenius norm o...
static const MagnitudeType blk_tol_default_
Block reorthogonalization threshold (default).
int max_ortho_steps_
Max number of (re)orthogonalization steps, including the first.
Class which defines basic traits for the operator type.
ScalarTraits< ScalarType >::magnitudeType normFrobenius() const
MagnitudeType sing_tol_
Singular block detection threshold.
Traits class which defines basic operations on multivectors.
static const MagnitudeType sing_tol_fast_
Singular block detection threshold (fast).
void innerProd(const MV &X, const MV &Y, Teuchos::SerialDenseMatrix< int, ScalarType > &Z) const
Provides the inner product defining the orthogonality concepts, using the provided operator...
bool blkOrtho1(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Routine to compute the block orthogonalization.
int findBasis(MV &X, Teuchos::RCP< MV > MX, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > C, bool completeBasis, int howMany=-1) const
Routine to find an orthonormal basis for X.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
bool is_null(const RCP< T > &p)
~ICGSOrthoManager()
Destructor.
Teuchos::ScalarTraits< ScalarType >::magnitudeType orthonormError(const MV &X) const
This method computes the error in orthonormality of a multivector, measured as the Frobenius norm of ...
MagnitudeType blk_tol_
Block reorthogonalization threshold.
Exception thrown to signal error in an orthogonalization manager method.
Teuchos::ScalarTraits< ScalarType > SCT
Teuchos::ScalarTraits< MagnitudeType > MGT
void setMyParamList(const RCP< ParameterList > ¶mList)
void resize(size_type new_size, const value_type &x=value_type())
TypeTo as(const TypeFrom &t)
RCP< ParameterList > getNonconstParameterList()
int normalize(MV &X, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B) const
This method calls normalize(X,Teuchos::null,B); see documentation for that function.
void setSingTol(const MagnitudeType sing_tol)
Set parameter for singular block detection.
MagnitudeType getSingTol() const
Return parameter for singular block detection.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using multiple steps ...
void push_back(const value_type &x)
Teuchos::RCP< Teuchos::ParameterList > defaultParams_
Default parameter list.
void project(MV &X, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
This method calls project(X,Teuchos::null,C,Q); see documentation for that function.
OperatorTraits< ScalarType, MV, OP > OPT
static const MagnitudeType blk_tol_fast_
Block reorthogonalization threshold (fast).
virtual int projectAndNormalizeWithMxImpl(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) 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
"Fast" but possibly unsafe or less accurate parameters.
bool blkOrtho(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Routine to compute the block orthogonalization.
Class which defines basic traits for the operator type.
Belos's templated virtual class for providing routines for orthogonalization and orthonormzalition of...
static const int max_ortho_steps_fast_
Max number of (re)orthogonalization steps, including the first (fast).
std::string label_
Label for timers.
Belos header file which uses auto-configuration information to include necessary C++ headers...
void project(MV &X, Teuchos::RCP< MV > MX, Teuchos::Array< Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > > C, Teuchos::ArrayView< Teuchos::RCP< const MV > > Q) const
Given a list of (mutually and internally) orthonormal bases Q, this method takes a multivector X and ...
int normalize(MV &X, Teuchos::RCP< MV > MX, Teuchos::RCP< Teuchos::SerialDenseMatrix< int, ScalarType > > B) const
This method takes a multivector X and attempts to compute an orthonormal basis for ...
Teuchos::RCP< Teuchos::ParameterList > getICGSFastParameters()
"Fast" but possibly unsafe or less accurate parameters.
Templated virtual class for providing orthogonalization/orthonormalization methods with matrix-based ...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const