MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
BelosXpetraStatusTestGenResSubNorm.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // MueLu: A package for multigrid based preconditioning
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 // Tobias Wiesner (tawiesn@sandia.gov)
43 //
44 // ***********************************************************************
45 //
46 // @HEADER
47 #ifndef BELOS_XPETRA_STATUS_TEST_GEN_RES_SUB_NORM_HPP
48 #define BELOS_XPETRA_STATUS_TEST_GEN_RES_SUB_NORM_HPP
49 
50 #include "Xpetra_ConfigDefs.hpp"
51 
53 
54 #include "MueLu_Exceptions.hpp"
55 
56 #include <BelosConfigDefs.hpp>
57 #include <BelosTypes.hpp>
58 #include <BelosOperatorT.hpp>
59 #include <BelosXpetraAdapterOperator.hpp>
60 #include <BelosStatusTestGenResSubNorm.hpp>
61 
62 
63 namespace Belos {
64 
68 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
69 class StatusTestGenResSubNorm<Scalar,Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>,Belos::OperatorT<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > >
70  : public StatusTestResNorm<Scalar,Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>,Belos::OperatorT<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > > {
71 
72  public:
73  // Convenience typedefs
77  typedef Belos::OperatorT<MV> OP;
78 
81  typedef MultiVecTraits<Scalar,MV> MVT;
82  typedef OperatorTraits<Scalar,MV,OP> OT;
83 
85 
86 
99  StatusTestGenResSubNorm( MagnitudeType Tolerance, size_t subIdx, int quorum = -1, bool showMaxResNormOnly = false )
100  : tolerance_(Tolerance),
101  subIdx_(subIdx),
102  quorum_(quorum),
103  showMaxResNormOnly_(showMaxResNormOnly),
104  resnormtype_(TwoNorm),
105  scaletype_(NormOfInitRes),
106  scalenormtype_(TwoNorm),
107  scalevalue_(Teuchos::ScalarTraits<MagnitudeType>::one ()),
108  status_(Undefined),
109  curBlksz_(0),
110  curNumRHS_(0),
111  curLSNum_(0),
112  numrhs_(0),
113  firstcallCheckStatus_(true),
114  firstcallDefineResForm_(true),
115  firstcallDefineScaleForm_(true),
116  mapExtractor_(Teuchos::null) { }
117 
119  virtual ~StatusTestGenResSubNorm() { };
121 
123 
124 
126 
132  int defineResForm(NormType TypeOfNorm) {
133  TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineResForm_==false,StatusTestError,
134  "StatusTestGenResSubNorm::defineResForm(): The residual form has already been defined.");
135  firstcallDefineResForm_ = false;
136 
137  resnormtype_ = TypeOfNorm;
138 
139  return(0);
140  }
141 
143 
163  int defineScaleForm( ScaleType TypeOfScaling, NormType TypeOfNorm, MagnitudeType ScaleValue = Teuchos::ScalarTraits<MagnitudeType>::one()) {
164  TEUCHOS_TEST_FOR_EXCEPTION(firstcallDefineScaleForm_==false,StatusTestError,
165  "StatusTestGenResSubNorm::defineScaleForm(): The scaling type has already been defined.");
166  firstcallDefineScaleForm_ = false;
167 
168  scaletype_ = TypeOfScaling;
169  scalenormtype_ = TypeOfNorm;
170  scalevalue_ = ScaleValue;
171 
172  return(0);
173  }
174 
176 
179  int setTolerance(MagnitudeType tolerance) {tolerance_ = tolerance; return(0);}
180 
182 
184  int setSubIdx ( size_t subIdx ) { subIdx_ = subIdx; return(0);}
185 
188  int setQuorum(int quorum) {quorum_ = quorum; return(0);}
189 
191  int setShowMaxResNormOnly(bool showMaxResNormOnly) {showMaxResNormOnly_ = showMaxResNormOnly; return(0);}
192 
194 
196 
197 
204  StatusType checkStatus(Iteration<Scalar,MV,OP>* iSolver) {
205  MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
206  const LinearProblem<Scalar,MV,OP>& lp = iSolver->getProblem();
207  // Compute scaling term (done once for each block that's being solved)
208  if (firstcallCheckStatus_) {
209  StatusType status = firstCallCheckStatusSetup(iSolver);
210  if(status==Failed) {
211  status_ = Failed;
212  return(status_);
213  }
214  }
215 
216  //
217  // This section computes the norm of the residual std::vector
218  //
219  if ( curLSNum_ != lp.getLSNumber() ) {
220  //
221  // We have moved on to the next rhs block
222  //
223  curLSNum_ = lp.getLSNumber();
224  curLSIdx_ = lp.getLSIndex();
225  curBlksz_ = (int)curLSIdx_.size();
226  int validLS = 0;
227  for (int i=0; i<curBlksz_; ++i) {
228  if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
229  validLS++;
230  }
231  curNumRHS_ = validLS;
232  curSoln_ = Teuchos::null;
233  //
234  } else {
235  //
236  // We are in the same rhs block, return if we are converged
237  //
238  if (status_==Passed) { return status_; }
239  }
240 
241  //
242  // Request the true residual for this block of right-hand sides.
243  //
244  Teuchos::RCP<MV> cur_update = iSolver->getCurrentUpdate();
245  curSoln_ = lp.updateSolution( cur_update );
246  Teuchos::RCP<MV> cur_res = MVT::Clone( *curSoln_, MVT::GetNumberVecs( *curSoln_ ) );
247  lp.computeCurrResVec( &*cur_res, &*curSoln_ );
248  std::vector<MagnitudeType> tmp_resvector( MVT::GetNumberVecs( *cur_res ) );
249  MvSubNorm( *cur_res, subIdx_, tmp_resvector, resnormtype_ );
250 
251  typename std::vector<int>::iterator p = curLSIdx_.begin();
252  for (int i=0; p<curLSIdx_.end(); ++p, ++i) {
253  // Check if this index is valid
254  if (*p != -1)
255  resvector_[*p] = tmp_resvector[i];
256  }
257 
258  //
259  // Compute the new linear system residuals for testing.
260  // (if any of them don't meet the tolerance or are NaN, then we exit with that status)
261  //
262  if ( scalevector_.size() > 0 ) {
263  typename std::vector<int>::iterator pp = curLSIdx_.begin();
264  for (; pp<curLSIdx_.end(); ++pp) {
265  // Check if this index is valid
266  if (*pp != -1) {
267  // Scale the std::vector accordingly
268  if ( scalevector_[ *pp ] != zero ) {
269  // Don't intentionally divide by zero.
270  testvector_[ *pp ] = resvector_[ *pp ] / scalevector_[ *pp ] / scalevalue_;
271  } else {
272  testvector_[ *pp ] = resvector_[ *pp ] / scalevalue_;
273  }
274  }
275  }
276  }
277  else {
278  typename std::vector<int>::iterator pp = curLSIdx_.begin();
279  for (; pp<curLSIdx_.end(); ++pp) {
280  // Check if this index is valid
281  if (*pp != -1)
282  testvector_[ *pp ] = resvector_[ *pp ] / scalevalue_;
283  }
284  }
285  // Check status of new linear system residuals and see if we have the quorum.
286  int have = 0;
287  ind_.resize( curLSIdx_.size() );
288  typename std::vector<int>::iterator p2 = curLSIdx_.begin();
289  for (; p2<curLSIdx_.end(); ++p2) {
290  // Check if this index is valid
291  if (*p != -1) {
292  // Check if any of the residuals are larger than the tolerance.
293  if (testvector_[ *p2 ] > tolerance_) {
294  // do nothing.
296  reset();
297  } else if (testvector_[ *p2 ] <= tolerance_) {
298  ind_[have] = *p2;
299  have++;
300  } else {
301  // Throw an std::exception if a NaN is found.
302  status_ = Failed;
303  TEUCHOS_TEST_FOR_EXCEPTION(true,StatusTestError,"StatusTestGenResSubNorm::checkStatus(): NaN has been detected.");
304  }
305  }
306  }
307  ind_.resize(have);
308  int need = (quorum_ == -1) ? curNumRHS_: quorum_;
309  status_ = (have >= need) ? Passed : Failed;
310  // Return the current status
311  return status_;
312  }
313 
315  StatusType getStatus() const {return(status_);};
317 
319 
320 
322  void reset() {
323  status_ = Undefined;
324  curBlksz_ = 0;
325  curLSNum_ = 0;
326  curLSIdx_.resize(0);
327  numrhs_ = 0;
328  ind_.resize(0);
329  firstcallCheckStatus_ = true;
330  curSoln_ = Teuchos::null;
331  }
332 
334 
336 
337 
339  void print(std::ostream& os, int indent = 0) const {
340  os.setf(std::ios_base::scientific);
341  for (int j = 0; j < indent; j ++)
342  os << ' ';
343  printStatus(os, status_);
344  os << resFormStr();
345  if (status_==Undefined)
346  os << ", tol = " << tolerance_ << std::endl;
347  else {
348  os << std::endl;
349  if(showMaxResNormOnly_ && curBlksz_ > 1) {
350  const MagnitudeType maxRelRes = *std::max_element(
351  testvector_.begin()+curLSIdx_[0],testvector_.begin()+curLSIdx_[curBlksz_-1]
352  );
353  for (int j = 0; j < indent + 13; j ++)
354  os << ' ';
355  os << "max{residual["<<curLSIdx_[0]<<"..."<<curLSIdx_[curBlksz_-1]<<"]} = " << maxRelRes
356  << ( maxRelRes <= tolerance_ ? " <= " : " > " ) << tolerance_ << std::endl;
357  }
358  else {
359  for ( int i=0; i<numrhs_; i++ ) {
360  for (int j = 0; j < indent + 13; j ++)
361  os << ' ';
362  os << "residual [ " << i << " ] = " << testvector_[ i ];
363  os << ((testvector_[i]<tolerance_) ? " < " : (testvector_[i]==tolerance_) ? " == " : (testvector_[i]>tolerance_) ? " > " : " " ) << tolerance_ << std::endl;
364  }
365  }
366  }
367  os << std::endl;
368  }
369 
371  void printStatus(std::ostream& os, StatusType type) const {
372  os << std::left << std::setw(13) << std::setfill('.');
373  switch (type) {
374  case Passed:
375  os << "Converged";
376  break;
377  case Failed:
378  os << "Unconverged";
379  break;
380  case Undefined:
381  default:
382  os << "**";
383  break;
384  }
385  os << std::left << std::setfill(' ');
386  return;
387  }
389 
391 
392 
394  Teuchos::RCP<MV> getSolution() { return curSoln_; }
395 
398  int getQuorum() const { return quorum_; }
399 
401  size_t getSubIdx() const { return subIdx_; }
402 
404  bool getShowMaxResNormOnly() { return showMaxResNormOnly_; }
405 
407  std::vector<int> convIndices() { return ind_; }
408 
410  MagnitudeType getTolerance() const {return(tolerance_);};
411 
413  const std::vector<MagnitudeType>* getTestValue() const {return(&testvector_);};
414 
416  const std::vector<MagnitudeType>* getResNormValue() const {return(&resvector_);};
417 
419  const std::vector<MagnitudeType>* getScaledNormValue() const {return(&scalevector_);};
420 
423  bool getLOADetected() const { return false; }
424 
426 
427 
430 
436  StatusType firstCallCheckStatusSetup(Iteration<Scalar,MV,OP>* iSolver) {
437  int i;
438  MagnitudeType zero = Teuchos::ScalarTraits<MagnitudeType>::zero();
439  MagnitudeType one = Teuchos::ScalarTraits<MagnitudeType>::one();
440  const LinearProblem<Scalar,MV,OP>& lp = iSolver->getProblem();
441  // Compute scaling term (done once for each block that's being solved)
442  if (firstcallCheckStatus_) {
443  //
444  // Get some current solver information.
445  //
446  firstcallCheckStatus_ = false;
447 
448  // try to access the underlying blocked operator
449  Teuchos::RCP<const OP> Op = lp.getOperator();
451  Teuchos::rcp_dynamic_cast<const Belos::XpetraOp<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
452  TEUCHOS_TEST_FOR_EXCEPTION(xOp.is_null(), MueLu::Exceptions::BadCast, "Bad cast from \'const Belos::OperatorT\' to \'const Belos::XpetraOp\'. The origin type is " << typeid(const OP).name() << ".");
454  xOp->getOperator();
455  TEUCHOS_TEST_FOR_EXCEPTION(xIntOp.is_null(), MueLu::Exceptions::BadCast, "Cannot access Xpetra::Operator stored in Belos::XpetraOperator.");
457  Teuchos::rcp_dynamic_cast<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(xIntOp);
458  TEUCHOS_TEST_FOR_EXCEPTION(xMat.is_null(), MueLu::Exceptions::RuntimeError, "Cannot access Xpetra::Matrix stored in Belos::XpetraOp. Error.");
460  TEUCHOS_TEST_FOR_EXCEPTION(bMat.is_null(), MueLu::Exceptions::BadCast, "Bad cast from \'const Xpetra::Matrix\' to \'const Xpetra::BlockedCrsMatrix\'. The origin type is " << typeid(const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>).name() << ". Note: you need a BlockedCrsMatrix object for the StatusTestGenResSubNorm to work!");
461  mapExtractor_ = bMat->getRangeMapExtractor();
462  TEUCHOS_TEST_FOR_EXCEPTION(mapExtractor_.is_null(), MueLu::Exceptions::RuntimeError, "Could not extract map extractor from BlockedCrsMatrix. Error.");
463  TEUCHOS_TEST_FOR_EXCEPTION(mapExtractor_->NumMaps()<=subIdx_, MueLu::Exceptions::RuntimeError, "The multivector is only split into " << mapExtractor_->NumMaps() << " sub parts. Cannot access sub-block " << subIdx_ << ".");
464 
465  // calculate initial norms
466  if (scaletype_== NormOfRHS) {
467  Teuchos::RCP<const MV> rhs = lp.getRHS();
468  numrhs_ = MVT::GetNumberVecs( *rhs );
469  scalevector_.resize( numrhs_ );
470  MvSubNorm( *rhs, subIdx_, scalevector_, scalenormtype_ );
471  }
472  else if (scaletype_==NormOfInitRes) {
473  Teuchos::RCP<const MV> init_res = lp.getInitResVec();
474  numrhs_ = MVT::GetNumberVecs( *init_res );
475  scalevector_.resize( numrhs_ );
476  MvSubNorm( *init_res, subIdx_, scalevector_, scalenormtype_ );
477  }
478  else if (scaletype_==NormOfPrecInitRes) {
479  Teuchos::RCP<const MV> init_res = lp.getInitPrecResVec();
480  numrhs_ = MVT::GetNumberVecs( *init_res );
481  scalevector_.resize( numrhs_ );
482  MvSubNorm( *init_res, subIdx_, scalevector_, scalenormtype_ );
483  }
484  else if (scaletype_==NormOfFullInitRes) {
485  Teuchos::RCP<const MV> init_res = lp.getInitResVec();
486  numrhs_ = MVT::GetNumberVecs( *init_res );
487  scalevector_.resize( numrhs_ );
488  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
489  scalevalue_ = one;
490  }
491  else if (scaletype_==NormOfFullPrecInitRes) {
492  Teuchos::RCP<const MV> init_res = lp.getInitPrecResVec();
493  numrhs_ = MVT::GetNumberVecs( *init_res );
494  scalevector_.resize( numrhs_ );
495  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
496  scalevalue_ = one;
497  }
498  else if (scaletype_==NormOfFullScaledInitRes) {
499  Teuchos::RCP<const MV> init_res = lp.getInitResVec();
500  numrhs_ = MVT::GetNumberVecs( *init_res );
501  scalevector_.resize( numrhs_ );
502  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
503  MvScalingRatio( *init_res, subIdx_, scalevalue_ );
504  }
505  else if (scaletype_==NormOfFullScaledPrecInitRes) {
506  Teuchos::RCP<const MV> init_res = lp.getInitPrecResVec();
507  numrhs_ = MVT::GetNumberVecs( *init_res );
508  scalevector_.resize( numrhs_ );
509  MVT::MvNorm( *init_res, scalevector_, scalenormtype_ );
510  MvScalingRatio( *init_res, subIdx_, scalevalue_ );
511  }
512  else {
513  numrhs_ = MVT::GetNumberVecs( *(lp.getRHS()) );
514  }
515 
516  resvector_.resize( numrhs_ );
517  testvector_.resize( numrhs_ );
518 
519  curLSNum_ = lp.getLSNumber();
520  curLSIdx_ = lp.getLSIndex();
521  curBlksz_ = (int)curLSIdx_.size();
522  int validLS = 0;
523  for (i=0; i<curBlksz_; ++i) {
524  if (curLSIdx_[i] > -1 && curLSIdx_[i] < numrhs_)
525  validLS++;
526  }
527  curNumRHS_ = validLS;
528  //
529  // Initialize the testvector.
530  for (i=0; i<numrhs_; i++) { testvector_[i] = one; }
531 
532  // Return an error if the scaling is zero.
533  if (scalevalue_ == zero) {
534  return Failed;
535  }
536  }
537  return Undefined;
538  }
540 
543 
545  std::string description() const
546  {
547  std::ostringstream oss;
548  oss << "Belos::StatusTestGenResSubNorm<>: " << resFormStr();
549  oss << ", tol = " << tolerance_;
550  return oss.str();
551  }
553 
554  protected:
555 
556  private:
557 
559 
560 
561  std::string resFormStr() const
562  {
563  std::ostringstream oss;
564  oss << "(";
565  oss << ((resnormtype_==OneNorm) ? "1-Norm" : (resnormtype_==TwoNorm) ? "2-Norm" : "Inf-Norm");
566  oss << " Exp";
567  oss << " Res Vec [" << subIdx_ << "]) ";
568 
569  // If there is no residual scaling, return current string.
570  if (scaletype_!=None)
571  {
572  // Insert division sign.
573  oss << "/ ";
574 
575  // Determine output string for scaling, if there is any.
576  if (scaletype_==UserProvided)
577  oss << " (User Scale)";
578  else {
579  oss << "(";
580  oss << ((scalenormtype_==OneNorm) ? "1-Norm" : (resnormtype_==TwoNorm) ? "2-Norm" : "Inf-Norm");
581  if (scaletype_==NormOfInitRes)
582  oss << " Res0 [" << subIdx_ << "]";
583  else if (scaletype_==NormOfPrecInitRes)
584  oss << " Prec Res0 [" << subIdx_ << "]";
585  else if (scaletype_==NormOfFullInitRes)
586  oss << " Full Res0 [" << subIdx_ << "]";
587  else if (scaletype_==NormOfFullPrecInitRes)
588  oss << " Full Prec Res0 [" << subIdx_ << "]";
589  else if (scaletype_==NormOfFullScaledInitRes)
590  oss << " scaled Full Res0 [" << subIdx_ << "]";
591  else if (scaletype_==NormOfFullScaledPrecInitRes)
592  oss << " scaled Full Prec Res0 [" << subIdx_ << "]";
593  else
594  oss << " RHS [" << subIdx_ << "]";
595  oss << ")";
596  }
597  }
598 
599  // TODO add a tagging name
600 
601  return oss.str();
602  }
603 
605 
607 
608 
609  // calculate norm of partial multivector
610  void MvSubNorm( const MV& mv, size_t block, std::vector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType>& normVec, NormType type = TwoNorm) {
611 
612  Teuchos::RCP<const MV> input = Teuchos::rcpFromRef(mv);
613 
614  Teuchos::RCP<const MV> SubVec = mapExtractor_->ExtractVector(input, block);
615  typedef MultiVecTraits<Scalar, MV> MVT;
616  MVT::MvNorm(*SubVec,normVec,type);
617  }
618 
619  // calculate ration of sub-vector length to full vector length (for scalevalue_)
620  void MvScalingRatio( const MV& mv, size_t block, MagnitudeType& lengthRatio) {
621  Teuchos::RCP<const MV> input = Teuchos::rcpFromRef(mv);
622 
623  Teuchos::RCP<const MV> SubVec = mapExtractor_->ExtractVector(input, block);
624 
625  lengthRatio = Teuchos::as<MagnitudeType>(SubVec->getGlobalLength()) / Teuchos::as<MagnitudeType>(input->getGlobalLength());
626  }
628 
630 
631 
633  MagnitudeType tolerance_;
634 
636  size_t subIdx_;
637 
639  int quorum_;
640 
642  bool showMaxResNormOnly_;
643 
645  NormType resnormtype_;
646 
648  ScaleType scaletype_;
649 
651  NormType scalenormtype_;
652 
654  MagnitudeType scalevalue_;
655 
657  std::vector<MagnitudeType> scalevector_;
658 
660  std::vector<MagnitudeType> resvector_;
661 
663  std::vector<MagnitudeType> testvector_;
664 
666  std::vector<int> ind_;
667 
669  Teuchos::RCP<MV> curSoln_;
670 
672  StatusType status_;
673 
675  int curBlksz_;
676 
678  int curNumRHS_;
679 
681  std::vector<int> curLSIdx_;
682 
684  int curLSNum_;
685 
687  int numrhs_;
688 
690  bool firstcallCheckStatus_;
691 
693  bool firstcallDefineResForm_;
694 
696  bool firstcallDefineScaleForm_;
697 
699  Teuchos::RCP<const ME> mapExtractor_;
701 
702 };
703 
704 } // namespace Belos
705 
706 #endif /* BELOS_XPETRA_STATUS_TEST_GEN_RES_SUB_NORM_HPP */
Exception indicating invalid cast attempted.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Exception throws to report errors in the internal logical of the program.
bool is_null() const