Anasazi  Version of the Day
AnasaziEpetraAdapter.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Anasazi: Block Eigensolvers Package
5 // Copyright 2004 Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 // @HEADER
41 
46 #ifndef ANASAZI_EPETRA_ADAPTER_HPP
47 #define ANASAZI_EPETRA_ADAPTER_HPP
48 
49 #include "AnasaziConfigDefs.hpp"
50 #include "Anasaziepetra_DLLExportMacro.h"
51 #include "AnasaziTypes.hpp"
52 #include "AnasaziMultiVec.hpp"
53 #include "AnasaziOperator.hpp"
55 
56 #include "Teuchos_Assert.hpp"
58 #include "Teuchos_FancyOStream.hpp"
59 
60 #include "Epetra_MultiVector.h"
61 #include "Epetra_Vector.h"
62 #include "Epetra_Operator.h"
63 #include "Epetra_Map.h"
64 #include "Epetra_LocalMap.h"
65 #include "Epetra_Comm.h"
66 
67 #if defined(HAVE_ANASAZI_TPETRA) && defined(HAVE_ANASAZI_TSQR)
68 # include <Tpetra_ConfigDefs.hpp> // HAVE_TPETRA_EPETRA
69 # if defined(HAVE_TPETRA_EPETRA)
70 # include <Epetra_TsqrAdaptor.hpp>
71 # endif // defined(HAVE_TPETRA_EPETRA)
72 #endif // defined(HAVE_ANASAZI_TPETRA) && defined(HAVE_ANASAZI_TSQR)
73 
74 namespace Anasazi {
75 
77 
78 
82  class EpetraMultiVecFailure : public AnasaziError {public:
83  EpetraMultiVecFailure(const std::string& what_arg) : AnasaziError(what_arg)
84  {}};
85 
89  class EpetraOpFailure : public AnasaziError {public:
90  EpetraOpFailure(const std::string& what_arg) : AnasaziError(what_arg)
91  {}};
92 
94 
96 
97 
102 
103  public:
106 
108  virtual Epetra_MultiVector* GetEpetraMultiVec() { return 0; }
109 
111  virtual const Epetra_MultiVector* GetEpetraMultiVec() const { return 0; }
112  };
113 
115 
117  //
118  //--------template class AnasaziEpetraMultiVec-----------------
119  //
121 
128  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraMultiVec : public MultiVec<double>, public Epetra_MultiVector, public EpetraMultiVecAccessor {
129  public:
131 
132 
134 
139  EpetraMultiVec(const Epetra_BlockMap& Map_in, const int numvecs);
140 
142  EpetraMultiVec(const Epetra_MultiVector & P_vec);
143 
145 
153  EpetraMultiVec(const Epetra_BlockMap& Map_in, double * array, const int numvecs, const int stride=0);
154 
156 
162  EpetraMultiVec(Epetra_DataAccess CV, const Epetra_MultiVector& P_vec, const std::vector<int>& index);
163 
165  virtual ~EpetraMultiVec() {};
166 
168 
170 
171 
176  MultiVec<double> * Clone ( const int numvecs ) const;
177 
183  MultiVec<double> * CloneCopy () const;
184 
192  MultiVec<double> * CloneCopy ( const std::vector<int>& index ) const;
193 
201  MultiVec<double> * CloneViewNonConst ( const std::vector<int>& index );
202 
210  const MultiVec<double> * CloneView ( const std::vector<int>& index ) const;
211 
213 
215  ptrdiff_t GetGlobalLength () const
216  {
217  if ( Map().GlobalIndicesLongLong() )
218  return static_cast<ptrdiff_t>( GlobalLength64() );
219  else
220  return static_cast<ptrdiff_t>( GlobalLength() );
221  }
222 
225  int GetNumberVecs () const { return NumVectors(); }
226 
228 
230 
231 
233  void MvTimesMatAddMv ( double alpha, const MultiVec<double>& A,
235  double beta );
236 
239  void MvAddMv ( double alpha, const MultiVec<double>& A,
240  double beta, const MultiVec<double>& B);
241 
244  void MvTransMv ( double alpha, const MultiVec<double>& A, Teuchos::SerialDenseMatrix<int,double>& B
245 #ifdef HAVE_ANASAZI_EXPERIMENTAL
246  , ConjType conj = Anasazi::CONJ
247 #endif
248  ) const;
249 
252  void MvDot ( const MultiVec<double>& A, std::vector<double> &b
253 #ifdef HAVE_ANASAZI_EXPERIMENTAL
254  , ConjType conj = Anasazi::CONJ
255 #endif
256  ) const;
257 
260  void MvScale ( double alpha ) {
261  TEUCHOS_TEST_FOR_EXCEPTION( this->Scale( alpha )!=0, EpetraMultiVecFailure,
262  "Anasazi::EpetraMultiVec::MvScale call to Epetra_MultiVector::Scale() returned a nonzero value.");
263  }
264 
267  void MvScale ( const std::vector<double>& alpha );
268 
270 
272 
276  void MvNorm ( std::vector<double> & normvec ) const {
277  if (((int)normvec.size() >= GetNumberVecs()) ) {
278  TEUCHOS_TEST_FOR_EXCEPTION( this->Norm2(&normvec[0])!=0, EpetraMultiVecFailure,
279  "Anasazi::EpetraMultiVec::MvNorm call to Epetra_MultiVector::Norm2() returned a nonzero value.");
280  }
281  };
283 
285 
286 
291  void SetBlock ( const MultiVec<double>& A, const std::vector<int>& index );
292 
295  void MvRandom() {
297  "Anasazi::EpetraMultiVec::MvRandom call to Epetra_MultiVector::Random() returned a nonzero value.");
298  };
299 
302  void MvInit ( double alpha ) {
303  TEUCHOS_TEST_FOR_EXCEPTION( this->PutScalar( alpha )!=0, EpetraMultiVecFailure,
304  "Anasazi::EpetraMultiVec::MvInit call to Epetra_MultiVector::PutScalar() returned a nonzero value.");
305  };
306 
308 
309 
311  Epetra_MultiVector* GetEpetraMultiVec() { return this; };
312 
314  const Epetra_MultiVector* GetEpetraMultiVec() const { return this; };
315 
317 
319 
321 
323  void MvPrint( std::ostream& os ) const { os << *this << std::endl; };
325 
326  private:
327  };
328  //-------------------------------------------------------------
329 
331  //
332  //--------template class AnasaziEpetraOp---------------------
333  //
335 
342  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraOp : public virtual Operator<double> {
343  public:
345 
346 
349 
351  ~EpetraOp();
353 
355 
356 
360  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
362 
363  private:
364 //use pragmas to disable some false-positive warnings for windows
365 // sharedlibs export
366 #ifdef _MSC_VER
367 #pragma warning(push)
368 #pragma warning(disable:4251)
369 #endif
371 #ifdef _MSC_VER
372 #pragma warning(pop)
373 #endif
374  };
375  //-------------------------------------------------------------
376 
378  //
379  //--------template class AnasaziEpetraGenOp--------------------
380  //
382 
396  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraGenOp : public virtual Operator<double>, public virtual Epetra_Operator {
397  public:
399 
404  bool isAInverse = true );
405 
407  ~EpetraGenOp();
408 
410 
412  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
413 
415 
417  int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const;
418 
420 
422  int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const;
423 
425  const char* Label() const { return "Epetra_Operator applying A^{-1}M"; };
426 
428  bool UseTranspose() const { return (false); };
429 
431  int SetUseTranspose(bool /*UseTranspose_in*/) { return 0; };
432 
434  bool HasNormInf() const { return (false); };
435 
437  double NormInf() const { return (-1.0); };
438 
440  const Epetra_Comm& Comm() const { return Epetra_AOp->Comm(); };
441 
443  const Epetra_Map& OperatorDomainMap() const { return Epetra_AOp->OperatorDomainMap(); };
444 
446  const Epetra_Map& OperatorRangeMap() const { return Epetra_AOp->OperatorRangeMap(); };
447 
448  private:
449  bool isAInverse;
450 
451 //use pragmas to disable some false-positive warnings for windows
452 // sharedlibs export
453 #ifdef _MSC_VER
454 #pragma warning(push)
455 #pragma warning(disable:4251)
456 #endif
459 #ifdef _MSC_VER
460 #pragma warning(pop)
461 #endif
462  };
463 
465  //
466  //--------template class AnasaziEpetraSymOp--------------------
467  //
469 
482  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraSymOp : public virtual Operator<double>, public virtual Epetra_Operator {
483  public:
485 
487  EpetraSymOp(const Teuchos::RCP<Epetra_Operator> &Op, bool isTrans = false );
488 
490  ~EpetraSymOp();
491 
493 
495  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
496 
498 
500  int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const;
501 
503 
506  int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const;
507 
509  const char* Label() const { return "Epetra_Operator applying A^TA or AA^T"; };
510 
512  bool UseTranspose() const { return (false); };
513 
515  int SetUseTranspose(bool /*UseTranspose_in*/) { return 0; };
516 
518  bool HasNormInf() const { return (false); };
519 
521  double NormInf() const { return (-1.0); };
522 
524  const Epetra_Comm& Comm() const { return Epetra_Op->Comm(); };
525 
527  const Epetra_Map& OperatorDomainMap() const { return Epetra_Op->OperatorDomainMap(); };
528 
530  const Epetra_Map& OperatorRangeMap() const { return Epetra_Op->OperatorRangeMap(); };
531 
532  private:
533 
534 //use pragmas to disable false-positive warnings in generating windows sharedlib exports
535 #ifdef _MSC_VER
536 #pragma warning(push)
537 #pragma warning(disable:4251)
538 #endif
540 #ifdef _MSC_VER
541 #pragma warning(pop)
542 #endif
543 
544  bool isTrans_;
545  };
546 
547 
549  //
550  //--------template class AnasaziEpetraSymMVOp---------------------
551  //
553 
566  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraSymMVOp : public virtual Operator<double> {
567  public:
569 
572  bool isTrans = false );
573 
576 
578 
580  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
581 
582  private:
583 
584 //use pragmas to disable some false-positive warnings for windows
585 // sharedlibs export
586 #ifdef _MSC_VER
587 #pragma warning(push)
588 #pragma warning(disable:4251)
589 #endif
591  Teuchos::RCP<const Epetra_Map> MV_localmap;
593 #ifdef _MSC_VER
594 #pragma warning(pop)
595 #endif
596 
597  bool isTrans_;
598  };
599 
601  //
602  //--------template class AnasaziEpetraWSymMVOp---------------------
603  //
605 
618  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraWSymMVOp : public virtual Operator<double> {
619  public:
622  const Teuchos::RCP<Epetra_Operator> &OP );
623 
626 
628 
630  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
631 
632  private:
633 //use pragmas to disable some false-positive warnings for windows
634 // sharedlibs export
635 #ifdef _MSC_VER
636 #pragma warning(push)
637 #pragma warning(disable:4251)
638 #endif
642  Teuchos::RCP<const Epetra_Map> MV_localmap;
644 #ifdef _MSC_VER
645 #pragma warning(pop)
646 #endif
647  };
648 
650  //
651  //--------template class AnasaziEpetraW2SymMVOp---------------------
652  //
654 
667  class ANASAZIEPETRA_LIB_DLL_EXPORT EpetraW2SymMVOp : public virtual Operator<double> {
668  public:
671  const Teuchos::RCP<Epetra_Operator> &OP );
672 
675 
677 
679  void Apply ( const MultiVec<double>& X, MultiVec<double>& Y ) const;
680 
681  private:
682 //use pragmas to disable some false-positive warnings for windows
683 // sharedlibs export
684 #ifdef _MSC_VER
685 #pragma warning(push)
686 #pragma warning(disable:4251)
687 #endif
691  Teuchos::RCP<const Epetra_Map> MV_localmap;
693 #ifdef _MSC_VER
694 #pragma warning(pop)
695 #endif
696  };
697 
698 
700  //
701  // Implementation of the Anasazi::MultiVecTraits for Epetra::MultiVector.
702  //
704 
715  template<>
716  class MultiVecTraits<double, Epetra_MultiVector>
717  {
718  public:
719 
721 
722 
728  Clone (const Epetra_MultiVector& mv, const int outNumVecs)
729  {
730  TEUCHOS_TEST_FOR_EXCEPTION(outNumVecs <= 0, std::invalid_argument,
731  "Belos::MultiVecTraits<double, Epetra_MultiVector>::"
732  "Clone(mv, outNumVecs = " << outNumVecs << "): "
733  "outNumVecs must be positive.");
734  // FIXME (mfh 13 Jan 2011) Anasazi currently lets Epetra fill in
735  // the entries of the returned multivector with zeros, but Belos
736  // does not. We retain this different behavior for now, but the
737  // two versions will need to be reconciled.
738  return Teuchos::rcp (new Epetra_MultiVector (mv.Map(), outNumVecs));
739  }
740 
746  CloneCopy (const Epetra_MultiVector& mv)
747  {
748  return Teuchos::rcp (new Epetra_MultiVector (mv));
749  }
750 
757  CloneCopy (const Epetra_MultiVector& mv, const std::vector<int>& index)
758  {
759  const int inNumVecs = GetNumberVecs (mv);
760  const int outNumVecs = index.size();
761 
762  // Simple, inexpensive tests of the index vector.
763  TEUCHOS_TEST_FOR_EXCEPTION(outNumVecs == 0, std::invalid_argument,
764  "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
765  "CloneCopy(mv, index = {}): At least one vector must be"
766  " cloned from mv.");
767  if (outNumVecs > inNumVecs)
768  {
769  std::ostringstream os;
770  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
771  "CloneCopy(mv, index = {";
772  for (int k = 0; k < outNumVecs - 1; ++k)
773  os << index[k] << ", ";
774  os << index[outNumVecs-1] << "}): There are " << outNumVecs
775  << " indices to copy, but only " << inNumVecs << " columns of mv.";
776  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
777  }
778 #ifdef TEUCHOS_DEBUG
779  // In debug mode, we perform more expensive tests of the index
780  // vector, to ensure all the elements are in range.
781  // Dereferencing the iterator is valid because index has length
782  // > 0.
783  const int minIndex = *std::min_element (index.begin(), index.end());
784  const int maxIndex = *std::max_element (index.begin(), index.end());
785 
786  if (minIndex < 0)
787  {
788  std::ostringstream os;
789  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
790  "CloneCopy(mv, index = {";
791  for (int k = 0; k < outNumVecs - 1; ++k)
792  os << index[k] << ", ";
793  os << index[outNumVecs-1] << "}): Indices must be nonnegative, but "
794  "the smallest index " << minIndex << " is negative.";
795  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
796  }
797  if (maxIndex >= inNumVecs)
798  {
799  std::ostringstream os;
800  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
801  "CloneCopy(mv, index = {";
802  for (int k = 0; k < outNumVecs - 1; ++k)
803  os << index[k] << ", ";
804  os << index[outNumVecs-1] << "}): Indices must be strictly less than "
805  "the number of vectors " << inNumVecs << " in mv; the largest index "
806  << maxIndex << " is out of bounds.";
807  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
808  }
809 #endif // TEUCHOS_DEBUG
810  // Cast to nonconst, because Epetra_MultiVector's constructor
811  // wants a nonconst int array argument. It doesn't actually
812  // change the entries of the array.
813  std::vector<int>& tmpind = const_cast< std::vector<int>& > (index);
814  return Teuchos::rcp (new Epetra_MultiVector (Epetra_DataAccess::Copy, mv, &tmpind[0], index.size()));
815  }
816 
818  CloneCopy (const Epetra_MultiVector& mv, const Teuchos::Range1D& index)
819  {
820  const int inNumVecs = GetNumberVecs (mv);
821  const int outNumVecs = index.size();
822  const bool validRange = outNumVecs > 0 && index.lbound() >= 0 &&
823  index.ubound() < inNumVecs;
824  if (! validRange)
825  {
826  std::ostringstream os;
827  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::Clone(mv,"
828  "index=[" << index.lbound() << ", " << index.ubound() << "]): ";
829  TEUCHOS_TEST_FOR_EXCEPTION(outNumVecs == 0, std::invalid_argument,
830  os.str() << "Column index range must be nonempty.");
831  TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
832  os.str() << "Column index range must be nonnegative.");
833  TEUCHOS_TEST_FOR_EXCEPTION(index.ubound() >= inNumVecs, std::invalid_argument,
834  os.str() << "Column index range must not exceed "
835  "number of vectors " << inNumVecs << " in the "
836  "input multivector.");
837  }
838  return Teuchos::rcp (new Epetra_MultiVector (Epetra_DataAccess::Copy, mv, index.lbound(), index.size()));
839  }
840 
847  CloneViewNonConst (Epetra_MultiVector& mv, const std::vector<int>& index)
848  {
849  const int inNumVecs = GetNumberVecs (mv);
850  const int outNumVecs = index.size();
851 
852  // Simple, inexpensive tests of the index vector.
853  TEUCHOS_TEST_FOR_EXCEPTION(outNumVecs == 0, std::invalid_argument,
854  "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
855  "CloneViewNonConst(mv, index = {}): The output view "
856  "must have at least one column.");
857  if (outNumVecs > inNumVecs)
858  {
859  std::ostringstream os;
860  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
861  "CloneViewNonConst(mv, index = {";
862  for (int k = 0; k < outNumVecs - 1; ++k)
863  os << index[k] << ", ";
864  os << index[outNumVecs-1] << "}): There are " << outNumVecs
865  << " indices to view, but only " << inNumVecs << " columns of mv.";
866  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
867  }
868 #ifdef TEUCHOS_DEBUG
869  // In debug mode, we perform more expensive tests of the index
870  // vector, to ensure all the elements are in range.
871  // Dereferencing the iterator is valid because index has length
872  // > 0.
873  const int minIndex = *std::min_element (index.begin(), index.end());
874  const int maxIndex = *std::max_element (index.begin(), index.end());
875 
876  if (minIndex < 0)
877  {
878  std::ostringstream os;
879  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
880  "CloneViewNonConst(mv, index = {";
881  for (int k = 0; k < outNumVecs - 1; ++k)
882  os << index[k] << ", ";
883  os << index[outNumVecs-1] << "}): Indices must be nonnegative, but "
884  "the smallest index " << minIndex << " is negative.";
885  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
886  }
887  if (maxIndex >= inNumVecs)
888  {
889  std::ostringstream os;
890  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::"
891  "CloneViewNonConst(mv, index = {";
892  for (int k = 0; k < outNumVecs - 1; ++k)
893  os << index[k] << ", ";
894  os << index[outNumVecs-1] << "}): Indices must be strictly less than "
895  "the number of vectors " << inNumVecs << " in mv; the largest index "
896  << maxIndex << " is out of bounds.";
897  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
898  }
899 #endif // TEUCHOS_DEBUG
900  // Cast to nonconst, because Epetra_MultiVector's constructor
901  // wants a nonconst int array argument. It doesn't actually
902  // change the entries of the array.
903  std::vector<int>& tmpind = const_cast< std::vector<int>& > (index);
904  return Teuchos::rcp (new Epetra_MultiVector (Epetra_DataAccess::View, mv, &tmpind[0], index.size()));
905  }
906 
908  CloneViewNonConst (Epetra_MultiVector& mv, const Teuchos::Range1D& index)
909  {
910  const bool validRange = index.size() > 0 &&
911  index.lbound() >= 0 &&
912  index.ubound() < mv.NumVectors();
913  if (! validRange)
914  {
915  std::ostringstream os;
916  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::CloneView"
917  "NonConst(mv,index=[" << index.lbound() << ", " << index.ubound()
918  << "]): ";
919  TEUCHOS_TEST_FOR_EXCEPTION(index.size() == 0, std::invalid_argument,
920  os.str() << "Column index range must be nonempty.");
921  TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
922  os.str() << "Column index range must be nonnegative.");
923  TEUCHOS_TEST_FOR_EXCEPTION(index.ubound() >= mv.NumVectors(),
924  std::invalid_argument,
925  os.str() << "Column index range must not exceed "
926  "number of vectors " << mv.NumVectors() << " in "
927  "the input multivector.");
928  }
929  return Teuchos::rcp (new Epetra_MultiVector (Epetra_DataAccess::View, mv, index.lbound(), index.size()));
930  }
931 
938  CloneView (const Epetra_MultiVector& mv, const std::vector<int>& index)
939  {
940  const int inNumVecs = GetNumberVecs (mv);
941  const int outNumVecs = index.size();
942 
943  // Simple, inexpensive tests of the index vector.
944  TEUCHOS_TEST_FOR_EXCEPTION(outNumVecs == 0, std::invalid_argument,
945  "Belos::MultiVecTraits<double,Epetra_MultiVector>::"
946  "CloneView(mv, index = {}): The output view "
947  "must have at least one column.");
948  if (outNumVecs > inNumVecs)
949  {
950  std::ostringstream os;
951  os << "Belos::MultiVecTraits<double,Epetra_MultiVector>::"
952  "CloneView(mv, index = {";
953  for (int k = 0; k < outNumVecs - 1; ++k)
954  os << index[k] << ", ";
955  os << index[outNumVecs-1] << "}): There are " << outNumVecs
956  << " indices to view, but only " << inNumVecs << " columns of mv.";
957  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
958  }
959 #ifdef TEUCHOS_DEBUG
960  // In debug mode, we perform more expensive tests of the index
961  // vector, to ensure all the elements are in range.
962  // Dereferencing the iterator is valid because index has length
963  // > 0.
964  const int minIndex = *std::min_element (index.begin(), index.end());
965  const int maxIndex = *std::max_element (index.begin(), index.end());
966 
967  if (minIndex < 0)
968  {
969  std::ostringstream os;
970  os << "Belos::MultiVecTraits<double,Epetra_MultiVector>::"
971  "CloneView(mv, index = {";
972  for (int k = 0; k < outNumVecs - 1; ++k)
973  os << index[k] << ", ";
974  os << index[outNumVecs-1] << "}): Indices must be nonnegative, but "
975  "the smallest index " << minIndex << " is negative.";
976  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
977  }
978  if (maxIndex >= inNumVecs)
979  {
980  std::ostringstream os;
981  os << "Belos::MultiVecTraits<double,Epetra_MultiVector>::"
982  "CloneView(mv, index = {";
983  for (int k = 0; k < outNumVecs - 1; ++k)
984  os << index[k] << ", ";
985  os << index[outNumVecs-1] << "}): Indices must be strictly less than "
986  "the number of vectors " << inNumVecs << " in mv; the largest index "
987  << maxIndex << " is out of bounds.";
988  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
989  }
990 #endif // TEUCHOS_DEBUG
991  // Cast to nonconst, because Epetra_MultiVector's constructor
992  // wants a nonconst int array argument. It doesn't actually
993  // change the entries of the array.
994  std::vector<int>& tmpind = const_cast< std::vector<int>& > (index);
995  return Teuchos::rcp (new Epetra_MultiVector (Epetra_DataAccess::View, mv, &tmpind[0], index.size()));
996  }
997 
999  CloneView (const Epetra_MultiVector& mv, const Teuchos::Range1D& index)
1000  {
1001  const bool validRange = index.size() > 0 &&
1002  index.lbound() >= 0 &&
1003  index.ubound() < mv.NumVectors();
1004  if (! validRange)
1005  {
1006  std::ostringstream os;
1007  os << "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::CloneView"
1008  "(mv,index=[" << index.lbound() << ", " << index.ubound()
1009  << "]): ";
1010  TEUCHOS_TEST_FOR_EXCEPTION(index.size() == 0, std::invalid_argument,
1011  os.str() << "Column index range must be nonempty.");
1012  TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
1013  os.str() << "Column index range must be nonnegative.");
1014  TEUCHOS_TEST_FOR_EXCEPTION(index.ubound() >= mv.NumVectors(),
1015  std::invalid_argument,
1016  os.str() << "Column index range must not exceed "
1017  "number of vectors " << mv.NumVectors() << " in "
1018  "the input multivector.");
1019  }
1020  return Teuchos::rcp (new Epetra_MultiVector(Epetra_DataAccess::View, mv, index.lbound(), index.size()));
1021  }
1022 
1024 
1026 
1027 
1029  static ptrdiff_t GetGlobalLength( const Epetra_MultiVector& mv )
1030  {
1031  if (mv.Map().GlobalIndicesLongLong())
1032  return static_cast<ptrdiff_t>( mv.GlobalLength64() );
1033  else
1034  return static_cast<ptrdiff_t>( mv.GlobalLength() );
1035  }
1036 
1038  static int GetNumberVecs( const Epetra_MultiVector& mv )
1039  { return mv.NumVectors(); }
1040 
1041  static bool HasConstantStride( const Epetra_MultiVector& mv )
1042  { return mv.ConstantStride(); }
1044 
1046 
1047 
1050  static void MvTimesMatAddMv( double alpha, const Epetra_MultiVector& A,
1052  double beta, Epetra_MultiVector& mv )
1053  {
1054  Epetra_LocalMap LocalMap(B.numRows(), 0, mv.Map().Comm());
1055  Epetra_MultiVector B_Pvec(Epetra_DataAccess::View, LocalMap, B.values(), B.stride(), B.numCols());
1056 
1057  TEUCHOS_TEST_FOR_EXCEPTION( mv.Multiply( 'N', 'N', alpha, A, B_Pvec, beta )!=0, EpetraMultiVecFailure,
1058  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvTimesMatAddMv call to Epetra_MultiVector::Multiply() returned a nonzero value.");
1059  }
1060 
1063  static void MvAddMv( double alpha, const Epetra_MultiVector& A, double beta, const Epetra_MultiVector& B, Epetra_MultiVector& mv )
1064  {
1065  // epetra mv.Update(alpha,A,beta,B,gamma) will check
1066  // alpha == 0.0
1067  // and
1068  // beta == 0.0
1069  // and based on this will call
1070  // mv.Update(beta,B,gamma)
1071  // or
1072  // mv.Update(alpha,A,gamma)
1073  //
1074  // mv.Update(alpha,A,gamma)
1075  // will then check for one of
1076  // gamma == 0
1077  // or
1078  // gamma == 1
1079  // or
1080  // alpha == 1
1081  // in that order. however, it will not look for the combination
1082  // alpha == 1 and gamma = 0
1083  // which is a common use case when we wish to assign
1084  // mv = A (in which case alpha == 1, beta == gamma == 0)
1085  // or
1086  // mv = B (in which case beta == 1, alpha == gamma == 0)
1087  //
1088  // therefore, we will check for these use cases ourselves
1089  if (beta == 0.0) {
1090  if (alpha == 1.0) {
1091  // assign
1092  mv = A;
1093  }
1094  else {
1095  // single update
1096  TEUCHOS_TEST_FOR_EXCEPTION( mv.Update( alpha, A, 0.0 )!=0, EpetraMultiVecFailure,
1097  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvAddMv call to Epetra_MultiVector::Update(alpha,A,0.0) returned a nonzero value.");
1098  }
1099  }
1100  else if (alpha == 0.0) {
1101  if (beta == 1.0) {
1102  // assign
1103  mv = B;
1104  }
1105  else {
1106  // single update
1107  TEUCHOS_TEST_FOR_EXCEPTION( mv.Update( beta, B, 0.0 )!=0, EpetraMultiVecFailure,
1108  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvAddMv call to Epetra_MultiVector::Update(beta,B,0.0) returned a nonzero value.");
1109  }
1110  }
1111  else {
1112  // double update
1113  TEUCHOS_TEST_FOR_EXCEPTION( mv.Update( alpha, A, beta, B, 0.0 )!=0, EpetraMultiVecFailure,
1114  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvAddMv call to Epetra_MultiVector::Update(alpha,A,beta,B,0.0) returned a nonzero value.");
1115  }
1116  }
1117 
1120  static void MvTransMv( double alpha, const Epetra_MultiVector& A, const Epetra_MultiVector& mv, Teuchos::SerialDenseMatrix<int,double>& B
1121 #ifdef HAVE_ANASAZI_EXPERIMENTAL
1122  , ConjType conj = Anasazi::CONJ
1123 #endif
1124  )
1125  {
1126  Epetra_LocalMap LocalMap(B.numRows(), 0, mv.Map().Comm());
1127  Epetra_MultiVector B_Pvec(Epetra_DataAccess::View, LocalMap, B.values(), B.stride(), B.numCols());
1128 
1129  TEUCHOS_TEST_FOR_EXCEPTION( B_Pvec.Multiply( 'T', 'N', alpha, A, mv, 0.0 )!=0, EpetraMultiVecFailure,
1130  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvTransMv call to Epetra_MultiVector::Multiply() returned a nonzero value.");
1131  }
1132 
1135  static void MvDot( const Epetra_MultiVector& A, const Epetra_MultiVector& B, std::vector<double> &b
1136 #ifdef HAVE_ANASAZI_EXPERIMENTAL
1137  , ConjType conj = Anasazi::CONJ
1138 #endif
1139  )
1140  {
1141 #ifdef TEUCHOS_DEBUG
1142  TEUCHOS_TEST_FOR_EXCEPTION(A.NumVectors() != B.NumVectors(),std::invalid_argument,
1143  "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::MvDot(A,B,b): A and B must have the same number of vectors.");
1144  TEUCHOS_TEST_FOR_EXCEPTION(b.size() != (unsigned int)A.NumVectors(),std::invalid_argument,
1145  "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::MvDot(A,B,b): b must have room for all dot products.");
1146 #endif
1147  TEUCHOS_TEST_FOR_EXCEPTION( A.Dot( B, &b[0] )!=0, EpetraMultiVecFailure,
1148  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvDot(A,B,b) call to Epetra_MultiVector::Dot() returned a nonzero value.");
1149  }
1150 
1152 
1154 
1158  static void MvNorm( const Epetra_MultiVector& mv, std::vector<double> &normvec )
1159  {
1160 #ifdef TEUCHOS_DEBUG
1161  TEUCHOS_TEST_FOR_EXCEPTION((unsigned int)mv.NumVectors() != normvec.size(),std::invalid_argument,
1162  "Anasazi::MultiVecTraits<double,Epetra_MultiVector>::MvNorm(mv,normvec): normvec must be the same size of mv.");
1163 #endif
1164  TEUCHOS_TEST_FOR_EXCEPTION( mv.Norm2(&normvec[0])!=0, EpetraMultiVecFailure,
1165  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvNorm call to Epetra_MultiVector::Norm2() returned a nonzero value.");
1166  }
1167 
1169 
1171 
1172 
1174  static void
1175  SetBlock (const Epetra_MultiVector& A,
1176  const std::vector<int>& index,
1177  Epetra_MultiVector& mv)
1178  {
1179  const int inNumVecs = GetNumberVecs (A);
1180  const int outNumVecs = index.size();
1181 
1182  // FIXME (mfh 13 Jan 2011) Belos allows A to have more columns
1183  // than index.size(), in which case we just take the first
1184  // index.size() columns of A. Anasazi requires that A have the
1185  // same number of columns as index.size(). Changing Anasazi's
1186  // behavior should not break existing Anasazi solvers, but the
1187  // tests need to be done.
1188  if (inNumVecs != outNumVecs)
1189  {
1190  std::ostringstream os;
1191  os << "Belos::MultiVecTraits<double,Epetra_MultiVector>::"
1192  "SetBlock(A, mv, index = {";
1193  if (outNumVecs > 0)
1194  {
1195  for (int k = 0; k < outNumVecs - 1; ++k)
1196  os << index[k] << ", ";
1197  os << index[outNumVecs-1];
1198  }
1199  os << "}): A has only " << inNumVecs << " columns, but there are "
1200  << outNumVecs << " indices in the index vector.";
1201  TEUCHOS_TEST_FOR_EXCEPTION(true, std::invalid_argument, os.str());
1202  }
1203  // Make a view of the columns of mv indicated by the index std::vector.
1205 
1206  // View of columns [0, outNumVecs-1] of the source multivector A.
1207  // If A has fewer columns than mv_view, then create a view of
1208  // the first outNumVecs columns of A.
1210  if (outNumVecs == inNumVecs)
1211  A_view = Teuchos::rcpFromRef (A); // Const, non-owning RCP
1212  else
1213  A_view = CloneView (A, Teuchos::Range1D(0, outNumVecs - 1));
1214 
1215  // Assignment calls Epetra_MultiVector::Assign(), which deeply
1216  // copies the data directly, ignoring the underlying
1217  // Epetra_Map(s). If A and mv don't have the same data
1218  // distribution (Epetra_Map), this may result in incorrect or
1219  // undefined behavior. Epetra_MultiVector::Update() also
1220  // ignores the Epetra_Maps, so we might as well just use the
1221  // (perhaps slightly cheaper) Assign() method via operator=().
1222  *mv_view = *A_view;
1223  }
1224 
1225  static void
1226  SetBlock (const Epetra_MultiVector& A,
1227  const Teuchos::Range1D& index,
1228  Epetra_MultiVector& mv)
1229  {
1230  const int numColsA = A.NumVectors();
1231  const int numColsMv = mv.NumVectors();
1232  // 'index' indexes into mv; it's the index set of the target.
1233  const bool validIndex = index.lbound() >= 0 && index.ubound() < numColsMv;
1234  // We can't take more columns out of A than A has.
1235  const bool validSource = index.size() <= numColsA;
1236 
1237  if (! validIndex || ! validSource)
1238  {
1239  std::ostringstream os;
1240  os << "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::SetBlock"
1241  "(A, index=[" << index.lbound() << ", " << index.ubound() << "], "
1242  "mv): ";
1243  TEUCHOS_TEST_FOR_EXCEPTION(index.lbound() < 0, std::invalid_argument,
1244  os.str() << "Range lower bound must be nonnegative.");
1245  TEUCHOS_TEST_FOR_EXCEPTION(index.ubound() >= numColsMv, std::invalid_argument,
1246  os.str() << "Range upper bound must be less than "
1247  "the number of columns " << numColsA << " in the "
1248  "'mv' output argument.");
1249  TEUCHOS_TEST_FOR_EXCEPTION(index.size() > numColsA, std::invalid_argument,
1250  os.str() << "Range must have no more elements than"
1251  " the number of columns " << numColsA << " in the "
1252  "'A' input argument.");
1253  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Should never get here!");
1254  }
1255 
1256  // View of columns [index.lbound(), index.ubound()] of the
1257  // target multivector mv. We avoid view creation overhead by
1258  // only creating a view if the index range is different than [0,
1259  // (# columns in mv) - 1].
1261  if (index.lbound() == 0 && index.ubound()+1 == numColsMv)
1262  mv_view = Teuchos::rcpFromRef (mv); // Non-const, non-owning RCP
1263  else
1264  mv_view = CloneViewNonConst (mv, index);
1265 
1266  // View of columns [0, index.size()-1] of the source multivector
1267  // A. If A has fewer columns than mv_view, then create a view
1268  // of the first index.size() columns of A.
1270  if (index.size() == numColsA)
1271  A_view = Teuchos::rcpFromRef (A); // Const, non-owning RCP
1272  else
1273  A_view = CloneView (A, Teuchos::Range1D(0, index.size()-1));
1274 
1275  // Assignment calls Epetra_MultiVector::Assign(), which deeply
1276  // copies the data directly, ignoring the underlying
1277  // Epetra_Map(s). If A and mv don't have the same data
1278  // distribution (Epetra_Map), this may result in incorrect or
1279  // undefined behavior. Epetra_MultiVector::Update() also
1280  // ignores the Epetra_Maps, so we might as well just use the
1281  // (perhaps slightly cheaper) Assign() method via operator=().
1282  *mv_view = *A_view;
1283  }
1284 
1285  static void
1286  Assign (const Epetra_MultiVector& A,
1287  Epetra_MultiVector& mv)
1288  {
1289  const int numColsA = GetNumberVecs (A);
1290  const int numColsMv = GetNumberVecs (mv);
1291  if (numColsA > numColsMv)
1292  {
1293  std::ostringstream os;
1294  os << "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::Assign"
1295  "(A, mv): ";
1296  TEUCHOS_TEST_FOR_EXCEPTION(numColsA > numColsMv, std::invalid_argument,
1297  os.str() << "Input multivector 'A' has "
1298  << numColsA << " columns, but output multivector "
1299  "'mv' has only " << numColsMv << " columns.");
1300  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error, "Should never get here!");
1301  }
1302  // View of the first [0, numColsA-1] columns of mv.
1304  if (numColsMv == numColsA)
1305  mv_view = Teuchos::rcpFromRef (mv); // Non-const, non-owning RCP
1306  else // numColsMv > numColsA
1307  mv_view = CloneView (mv, Teuchos::Range1D(0, numColsA - 1));
1308 
1309  // Assignment calls Epetra_MultiVector::Assign(), which deeply
1310  // copies the data directly, ignoring the underlying
1311  // Epetra_Map(s). If A and mv don't have the same data
1312  // distribution (Epetra_Map), this may result in incorrect or
1313  // undefined behavior. Epetra_MultiVector::Update() also
1314  // ignores the Epetra_Maps, so we might as well just use the
1315  // (perhaps slightly cheaper) Assign() method via operator=().
1316  *mv_view = A;
1317  }
1318 
1321  static void MvScale ( Epetra_MultiVector& mv, double alpha )
1322  {
1323  TEUCHOS_TEST_FOR_EXCEPTION( mv.Scale( alpha )!=0, EpetraMultiVecFailure,
1324  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvScale call to Epetra_MultiVector::Scale(mv,double alpha) returned a nonzero value.");
1325  }
1326 
1329  static void MvScale ( Epetra_MultiVector& mv, const std::vector<double>& alpha )
1330  {
1331  // Check to make sure the vector is as long as the multivector has columns.
1332  int numvecs = mv.NumVectors();
1333 #ifdef TEUCHOS_DEBUG
1334  TEUCHOS_TEST_FOR_EXCEPTION( alpha.size() != (unsigned int)numvecs, std::invalid_argument,
1335  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvScale(mv,vector alpha): size of alpha inconsistent with number of vectors in mv.")
1336 #endif
1337  for (int i=0; i<numvecs; i++) {
1338  TEUCHOS_TEST_FOR_EXCEPTION( mv(i)->Scale(alpha[i])!=0, EpetraMultiVecFailure,
1339  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvScale call to Epetra_MultiVector::Scale() returned a nonzero value.");
1340  }
1341  }
1342 
1345  static void MvRandom( Epetra_MultiVector& mv )
1346  {
1348  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvRandom call to Epetra_MultiVector::Random() returned a nonzero value.");
1349  }
1350 
1353  static void MvInit( Epetra_MultiVector& mv, double alpha = Teuchos::ScalarTraits<double>::zero() )
1354  {
1355  TEUCHOS_TEST_FOR_EXCEPTION( mv.PutScalar(alpha)!=0, EpetraMultiVecFailure,
1356  "Anasazi::MultiVecTraits<double, Epetra_MultiVector>::MvInit call to Epetra_MultiVector::PutScalar() returned a nonzero value.");
1357  }
1358 
1360 
1362 
1363 
1366  static void MvPrint( const Epetra_MultiVector& mv, std::ostream& os )
1367  { os << mv << std::endl; }
1368 
1370 
1371 #if defined(HAVE_ANASAZI_TPETRA) && defined(HAVE_ANASAZI_TSQR)
1372 # if defined(HAVE_TPETRA_EPETRA)
1373  typedef Epetra::TsqrAdaptor tsqr_adaptor_type;
1379 # endif // defined(HAVE_TPETRA_EPETRA)
1380 #endif // defined(HAVE_ANASAZI_TPETRA) && defined(HAVE_ANASAZI_TSQR)
1381  };
1382 
1384  //
1385  // Implementation of the Anasazi::OperatorTraits for Epetra::Operator.
1386  //
1388 
1400  template <>
1401  class OperatorTraits < double, Epetra_MultiVector, Epetra_Operator >
1402  {
1403  public:
1404 
1408  static void Apply ( const Epetra_Operator& Op,
1409  const Epetra_MultiVector& x,
1410  Epetra_MultiVector& y )
1411  {
1412 #ifdef TEUCHOS_DEBUG
1413  TEUCHOS_TEST_FOR_EXCEPTION(x.NumVectors() != y.NumVectors(),std::invalid_argument,
1414  "Anasazi::OperatorTraits<double,Epetra_MultiVector,Epetra_Operator>::Apply(Op,x,y): x and y must have the same number of columns.");
1415 #endif
1416  int ret = Op.Apply(x,y);
1418  "Anasazi::OperatorTraits<double,Epetra_Multivector,Epetra_Operator>::Apply(): Error in Epetra_Operator::Apply(). Code " << ret);
1419  }
1420 
1421  };
1422 
1423  template<>
1424  struct OutputStreamTraits<Epetra_Operator>
1425  {
1427  getOutputStream (const Epetra_Operator& op, int rootRank = 0)
1428  {
1429  Teuchos::RCP<Teuchos::FancyOStream> fos = Teuchos::getFancyOStream(Teuchos::rcpFromRef(std::cout));
1430  const Epetra_Comm & comm = op.Comm();
1431 
1432  // Select minimum MPI rank as the root rank for printing if provided rank is less than 0.
1433  int myRank = comm.MyPID();
1434  int numProcs = comm.NumProc();
1435  if (rootRank < 0)
1436  {
1437  comm.MinAll( &myRank, &rootRank, 1 );
1438  }
1439 
1440  // This is irreversible, but that's only a problem if the input std::ostream
1441  // is actually a Teuchos::FancyOStream on which this method has been
1442  // called before, with a different root rank.
1443  fos->setProcRankAndSize (myRank, numProcs);
1444  fos->setOutputToRootOnly (rootRank);
1445  return fos;
1446  }
1447  };
1448 
1449 } // end of Anasazi namespace
1450 
1451 #endif
1452 // end of file ANASAZI_EPETRA_ADAPTER_HPP
basic_FancyOStream & setProcRankAndSize(const int procRank, const int numProcs)
void MvRandom()
Fill the vectors in *this with random numbers.
Adapter class for creating an operators often used in solving generalized eigenproblems.
EpetraMultiVecAccessor is an interfaceto allow any Anasazi::MultiVec implementation that is based on ...
const Epetra_Comm & Comm() const
Returns the Epetra_Comm communicator associated with this operator.
ScalarType * values() const
static void MvAddMv(double alpha, const Epetra_MultiVector &A, double beta, const Epetra_MultiVector &B, Epetra_MultiVector &mv)
Replace mv with .
void MvInit(double alpha)
Replace each element of the vectors in *this with alpha.
static void SetBlock(const Epetra_MultiVector &A, const std::vector< int > &index, Epetra_MultiVector &mv)
Copy the vectors in A to a set of vectors in mv indicated by the indices given in index...
virtual ~EpetraMultiVecAccessor()
Destructor.
Adapter class for creating a weighted symmetric operator from an Epetra_MultiVector and Epetra_Operat...
ptrdiff_t GetGlobalLength() const
The number of rows in the multivector.
static void MvInit(Epetra_MultiVector &mv, double alpha=Teuchos::ScalarTraits< double >::zero())
Replace each element of the vectors in mv with alpha.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Virtual base class which defines basic traits for the operator type.
const char * Label() const
Returns a character string describing the operator.
virtual Epetra_MultiVector * GetEpetraMultiVec()
Return the pointer to the Epetra_MultiVector object.
static void MvPrint(const Epetra_MultiVector &mv, std::ostream &os)
Print the mv multi-vector to the os output stream.
static void MvRandom(Epetra_MultiVector &mv)
Replace the vectors in mv with random vectors.
static void MvTransMv(double alpha, const Epetra_MultiVector &A, const Epetra_MultiVector &mv, Teuchos::SerialDenseMatrix< int, double > &B)
Compute a dense matrix B through the matrix-matrix multiply .
static void Assign(const MV &A, MV &mv)
mv := A
bool UseTranspose() const
Returns the current UseTranspose setting [always false for this operator].
An exception class parent to all Anasazi exceptions.
double NormInf() const
Returns the infinity norm of the global matrix [not functional for this operator].
Interface for multivectors used by Anasazi&#39; linear solvers.
Namespace Anasazi contains the classes, structs, enums and utilities used by the Anasazi package...
Basic adapter class for Anasazi::Operator that uses Epetra_Operator.
static void MvDot(const Epetra_MultiVector &A, const Epetra_MultiVector &B, std::vector< double > &b)
Compute a vector b where the components are the individual dot-products of the i-th columns of A and ...
bool HasNormInf() const
Returns true if this object can provide an approximate inf-norm [always false for this operator]...
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
double NormInf() const
Returns the infinity norm of the global matrix [not functional for this operator].
const Epetra_Comm & Comm() const
Returns the Epetra_Comm communicator associated with this operator.
ConjType
Enumerated types used to specify conjugation arguments.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
OrdinalType numRows() const
bool UseTranspose() const
Returns the current UseTranspose setting [always false for this operator].
static Teuchos::RCP< Epetra_MultiVector > CloneCopy(const Epetra_MultiVector &mv, const std::vector< int > &index)
Creates a new Epetra_MultiVector and copies the selected contents of mv into the new vector (deep cop...
Epetra_MultiVector * GetEpetraMultiVec()
Return the pointer to the Epetra_MultiVector object.
int GetNumberVecs() const
The number of vectors (i.e., columns) in the multivector.
const char * Label() const
Returns a character string describing the operator.
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
static void MvScale(Epetra_MultiVector &mv, const std::vector< double > &alpha)
Scale each element of the i-th vector in mv with alpha[i].
static void MvTimesMatAddMv(double alpha, const Epetra_MultiVector &A, const Teuchos::SerialDenseMatrix< int, double > &B, double beta, Epetra_MultiVector &mv)
Update mv with .
Traits class which defines basic operations on multivectors.
static Teuchos::RCP< MV > CloneCopy(const MV &mv)
Creates a new MV and copies contents of mv into the new vector (deep copy).
static Teuchos::RCP< const Epetra_MultiVector > CloneView(const Epetra_MultiVector &mv, const std::vector< int > &index)
Creates a new const Epetra_MultiVector that shares the selected contents of mv (shallow copy)...
Adapter class for creating a weighted operator from an Epetra_MultiVector and Epetra_Operator.
virtual const Epetra_MultiVector * GetEpetraMultiVec() const
Return the pointer to the Epetra_MultiVector object.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
Ordinal ubound() const
static Teuchos::RCP< Epetra_MultiVector > Clone(const Epetra_MultiVector &mv, const int outNumVecs)
Creates a new empty Epetra_MultiVector containing numVecs columns.
Anasazi header file which uses auto-configuration information to include necessary C++ headers...
basic_FancyOStream & setOutputToRootOnly(const int rootRank)
void MvPrint(std::ostream &os) const
Print *this EpetraMultiVec.
static Teuchos::RCP< MV > CloneViewNonConst(MV &mv, const std::vector< int > &index)
Creates a new MV that shares the selected contents of mv (shallow copy).
Ordinal lbound() const
static ptrdiff_t GetGlobalLength(const Epetra_MultiVector &mv)
Obtain the vector length of mv.
int SetUseTranspose(bool)
If set true, the transpose of this operator will be applied [not functional for this operator]...
Adapter class for creating a symmetric operator from an Epetra_Operator.
static Teuchos::RCP< const MV > CloneView(const MV &mv, const std::vector< int > &index)
Creates a new const MV that shares the selected contents of mv (shallow copy).
const Epetra_MultiVector * GetEpetraMultiVec() const
Return the pointer to the Epetra_MultiVector object.
void MvScale(double alpha)
Scale each element of the vectors in *this with alpha.
static void MvScale(Epetra_MultiVector &mv, double alpha)
Scale each element of the vectors in mv with alpha.
static void SetBlock(const MV &A, const std::vector< int > &index, MV &mv)
Copy the vectors in A to a set of vectors in mv indicated by the indices given in index...
static void Apply(const Epetra_Operator &Op, const Epetra_MultiVector &x, Epetra_MultiVector &y)
This method takes the Epetra_MultiVector x and applies the Epetra_Operator Op to it resulting in the ...
static int GetNumberVecs(const MV &mv)
Obtain the number of vectors in mv.
Templated virtual class for creating operators that can interface with the Anasazi::OperatorTraits cl...
static Teuchos::RCP< Epetra_MultiVector > CloneViewNonConst(Epetra_MultiVector &mv, const std::vector< int > &index)
Creates a new Epetra_MultiVector that shares the selected contents of mv (shallow copy)...
Ordinal size() const
void MvNorm(std::vector< double > &normvec) const
Compute the 2-norm of each individual vector of *this. Upon return, normvec[i] holds the 2-norm of t...
EpetraMultiVecFailure is thrown when a return value from an Epetra call on an Epetra_MultiVector is n...
Types and exceptions used within Anasazi solvers and interfaces.
Abstract class definition for Anasazi output stream.
OrdinalType stride() const
OrdinalType numCols() const
Interface for multivectors used by Anasazi&#39;s linear solvers.
Adapter class for creating a symmetric operator from an Epetra_MultiVector.
static void MvNorm(const Epetra_MultiVector &mv, std::vector< double > &normvec)
Compute the 2-norm of each individual vector of mv. Upon return, normvec[i] holds the value of ...
static int GetNumberVecs(const Epetra_MultiVector &mv)
Obtain the number of vectors in mv.
static Teuchos::RCP< Epetra_MultiVector > CloneCopy(const Epetra_MultiVector &mv)
Creates a new Epetra_MultiVector and copies contents of mv into the new vector (deep copy)...
Anasazi&#39;s templated virtual class for constructing an operator that can interface with the OperatorTr...
virtual ~EpetraMultiVec()
Destructor.
int SetUseTranspose(bool)
If set true, the transpose of this operator will be applied [not functional for this operator]...
EpetraOpFailure is thrown when a return value from an Epetra call on an Epetra_Operator is non-zero...
Exceptions thrown to signal error in operator application.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
Basic adapter class for Anasazi::MultiVec that uses Epetra_MultiVector.
bool HasNormInf() const
Returns true if this object can provide an approximate inf-norm [always false for this operator]...